      SUBROUTINE RWQC1
C
C CHANGE RECORD
C   Merged SNL and DSINTL codes
C READ IN FROM THE UNIT #8
C: I/O CONTROL VARIABLES
C: SPATIALLY AND TEMPORALLY CONSTANT REAL PARAMETERS
C
      USE GLOBAL
C
 	IMPLICIT NONE

      CHARACTER TITLE(5)*79, CCMRM*1
      CHARACTER LINE*255
      REAL,SAVE,ALLOCATABLE,DIMENSION(:)::XDSL
      REAL,SAVE,ALLOCATABLE,DIMENSION(:)::XPSL
      REAL      WQTDTEMP(1000),CONV1,CONV2,WQDIUDT,XC,XP,XPC,XPD,XPG
      REAL      XN,XNC,XND,XNG,WTEMP,TT20,WQTT,O2WQ_,WQTAMD,TVARWQ
      REAL      TWQTSB,TWQTSE,WQTSDT,XWQCHL,XWQTAMP,XWQPO4D,XWQSAD  ! C06
      REAL      WQTMC1,WQTMC2,WQTMD1,WQTMD2,WQTMG1,WQTMG2,WQTMP1,   ! C11
     &          WQTMP2
      REAL      WQKG1C,WQKG2C,WQKG1D,WQKG2D,WQKG1G,WQKG2G,WQKG1P,   ! C12
     &          WQKG2P
      REAL      WQTRC,WQTRD,WQTRM,WQKTBC,WQKTBD,WQKTBG              ! C13
      REAL      WQTRHDR,WQTRMNL,WQKTHDR,WQKTMNL                     ! C17
      REAL      WQTNIT,WQKN1,WQKN2                                  ! C25
      REAL      WQKSU,WQTRSUA,WQKTSUA                               ! C27
      REAL      WQTRCOD,WQKTCOD                                     ! C28
      REAL      WQBFTAM,WQTTAM,WQKFCB,WQTFCB,WQKTAM                 ! C29
      REAL      XMRM1, XMRM2, XMRM3,XMRM4,XMRMA,XMRMB,XMRMC,XMRMD,  ! MACROALGAE
     &          XMRME
      REAL      XPSQ,XPSQ1,XDSQ,XMUD
      REAL      FRACLAY,FRACLAYKM1,FHLAYC,FHLAYCKM1,WQVMTMP
      INTEGER   M,N1,II,JJ,KK,NT,ISSKIP,NW,ND,LF,LL,L
      INTEGER   IWQDT,IWQKIN,IZ,IN,IJKC,IWQZX,IZMUD,IZSAND
      INTEGER   I,J,K,MM
      INTEGER   ITMP !From C48 to tel where to add mass to irrigation
      INTEGER   I_TEMP, J_TEMP
      LOGICAL   CELL_INSIDE_DOMAIN
C
      PARAMETER (CONV1=1.0,CONV2=8.64E4)
C
      IF(.NOT.ALLOCATED(XDSL))THEN
		ALLOCATE(XDSL(NWQVM))
		ALLOCATE(XPSL(NWQVM))
	    XDSL=0.0
	    XPSL=0.0
      ENDIF
C
      OPEN(2,FILE='WQ3D'//ans(partid2)//'.OUT',STATUS='UNKNOWN',POSITION='APPEND')
      IF(PARTID==MASTER_TASK)PRINT *,'WQ: READING WQ3DWC.INP - MAIN WATER QUALITY CONTROL FILE'
      OPEN(1,FILE='WQ3DWC.INP',STATUS='UNKNOWN')
C
C READ FIRST LINE IN WQ3DWC.INP FILE.  IF FIRST CHARACTER IS '#', THEN
C THIS IS THE NEW VERSION WITH ANNOTATED COMMENTS ADDED (I.E., USES THE
C SKIPCOMM SUBROUTINE TO SKIP COMMENT LINES.  COMMENT LINES BEGIN WITH
C A "C", "C", OR "#" CHARACTER IN COLUMN 1.  IF "#" IS NOT FOUND AS THE
C FIRST CHARACTER IN THE FILE, THEN THE OLD METHOD OF READING THE
C WQ3DWC.INP FILE IS USED TO PRESERVE BACKWARD COMPATABILITY.
C
      ISSKIP = 0
      READ(1,'(A1)') CCMRM
      BACKSPACE(1)
      IF(CCMRM .EQ. '#') ISSKIP = 1
      CCMRM = '#'
C *** C01
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      READ(1,90) (TITLE(M), M=1,3)
      WRITE(2,90) (TITLE(M), M=1,3)
C
C I/O CONTROL VARIABLES
C      READ(1,999)
C
      WRITE(2,999)
C *** C02
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      READ(1,90) TITLE(1)
      WRITE(2,90) TITLE(1)
      WRITE(2,999)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      READ(1,*) ISWQLVL,NWQV,NWQZ,NWQPS_GL,NWQTD,NWQTS,NTSWQV,NSMG,NSMZ,NTSSMV,NSMTS,NWQKDPT
      WRITE(2,*) ISWQLVL,NWQV,NWQZ,NWQPS_GL,NWQTD,NWQTS,NTSWQV,NSMG,NSMZ,NTSSMV,NSMTS,NWQKDPT
      IF(ISWQLVL.LT.0.OR.ISWQLVL.GT.4)STOP 'BAD KINETICS OPTION'   ! *** PMC
C *** C02A
      ! *** ONLY USED WHEN ISWQLVL=1-3
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      READ(1,*)  (ISTRWQ(NW),NW=1,NWQV)
      WRITE(2,*) (ISTRWQ(NW),NW=1,NWQV)
      IF(ISWQLVL.EQ.0)THEN
        ISTRWQ(1:NWQV)=0
        ISTRWQ(6)=1
        ISTRWQ(14)=1
        ISTRWQ(19)=1
      ENDIF
C *** C03
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      READ(1,*)  IWQDT,IWQM,IWQBEN,IWQSI,IWQFCB,IWQSRP,IWQSTOX,IWQKA(1),IWQVLIM
      WRITE(2,*) IWQDT,IWQM,IWQBEN,IWQSI,IWQFCB,IWQSRP,IWQSTOX,IWQKA(1),IWQVLIM
C *** C04
      IF(ISSKIP .GT. 0)CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0)READ(1,999)
      READ(1,*)  IWQZ,IWQNC,IWQRST,NDMWQ,LDMWQ,NDDOAVG,NDLTAVG,IDNOTRVA
      WRITE(2,*) IWQZ,IWQNC,IWQRST,NDMWQ,LDMWQ,NDDOAVG,NDLTAVG,IDNOTRVA
C *** C05
      IF(ISSKIP .GT. 0)CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0)READ(1,999)
      READ(1,*)  IWQICI,IWQAGR,IWQSTL,IWQSUN,IWQPSL,IWQNPL,ISDIURDO,WQDIUDT,IWQKIN
      WRITE(2,*) IWQICI,IWQAGR,IWQSTL,IWQSUN,IWQPSL,IWQNPL,ISDIURDO,WQDIUDT,IWQKIN
      IWQDIUDT = NINT(WQDIUDT*3600.0/DT)
      WRITE(2,83)': FREQUENCY OF DIURNAL DO OUTPUT (IN DT UNIT) =',IWQDIUDT
      WRITE(2,83)'* IWQDT (DTWQ(D) = DT(S)*IWQDT/86400)        = ',IWQDT
      DTD = DT/86400.0
C      DTWQ = DTD*REAL(IWQDT)*REAL(NWQKDPT)  PMC
      DTWQ = DTD*REAL(NWQKDPT)
      DTWQO2 = DTWQ*0.5
      !IF(IWQM.EQ.1)THEN
      WRITE(2,80)'* FULL VERSION WITH 21 VARIABLES IS ACTIVATED     '
      !ELSE IF(IWQM.EQ.2)THEN
      !STOP 'SMALL VERSION WITH 9 VARIABLES IS NOT OPERATIONAL, STOPPING'
      !ELSE
      !STOP '** ERROR!!! INVALID IWQM VALUE **'
      !ENDIF
      IF(IWQBEN.EQ.1)THEN
        WRITE(2,80)'* SEDIMENT PROCESS MODEL IS ACTIVATED             '
      ELSEIF(IWQBEN.EQ.0)THEN
        WRITE(2,80)'* SPATIALLY/TEMPORALLY CONSTANT BF IS SPECIFIED   '
      ELSEIF(IWQBEN.EQ.2)THEN
        WRITE(2,80)'* SPATIALLY AND/OR TEMPORALLY-VARYING BF SPECIFIED'
      ELSE
        STOP '** ERROR!!! INVALID IWQBEN VALUE **'
      ENDIF
      IF(IWQSI.EQ.1)THEN
        WRITE(2,80)'* SILICA STATE VARIABLES (SU & SA) ARE MODELED    '
      ELSE
        WRITE(2,80)'* NO SILICA (SU & SU) LIMITATION                  '
      ENDIF
      IF(IWQFCB.EQ.1)THEN
        WRITE(2,80)'* FCB (FECAL COLIFORM BACTERIA) IS MODELED        '
      ELSE
        WRITE(2,80)'* FCB (FECAL COLIFORM BACTERIA) IS NOT MODELED    '
      ENDIF
      IF(IWQSRP.EQ.1)THEN
        WRITE(2,80)'* TAM IS USED FOR SORPTION OF PO4T/SA: MODEL TAM  '
      ELSEIF(IWQSRP.EQ.2)THEN
        WRITE(2,80)'* TSS IS USED FOR SORPTION OF PO4T/SA: MODEL TSS  '
      IF(NSED.EQ.0) STOP 'ERROR! INCOMPATIBLE ISTRAN(6)/NSED/IWQSRP'
!      IF(ISTRAN(6).NE.1) STOP 'ERROR! INCOMPATIBLE ISTRAN(6)/IWQSRP'  !SCJ changed to work with VEFDC for sorption to TSS
      ELSE
        WRITE(2,80)'* NO SORPTION OF PO4T/SA: MAY MODEL TSS & NO TAM  '
      ENDIF
      IF(IWQSTOX.EQ.1)THEN
        WRITE(2,80)'* SALINITY TOXICITY IS APPLIED TO CYANOBACTERIA   '
      ELSE
        WRITE(2,80)'* NO SALINITY TOXICITY: SALTWATER CYANOBACTERIA   '
      ENDIF
      IF(IWQKA(1).EQ.0)THEN
        WRITE(2,80)'* USER-SPECIFIED CONSTANT REAERATION SET TO WQKRO '
        WRITE(2,80)'*   REAERATION DUE TO WIND SET TO ZERO            '
      ELSEIF(IWQKA(1).EQ.1)THEN
        WRITE(2,80)'* USER-SPECIFIED CONSTANT REAERATION SET TO WQKRO '
        WRITE(2,80)'*   REAERATION DUE TO WIND ADDED TO WQKRO         '
      ELSEIF(IWQKA(1).EQ.2)THEN
        WRITE(2,80)'* OCONNOR-DOBBINS REAERATION FORMULA IS USED      '
      ELSEIF(IWQKA(1).EQ.3)THEN
        WRITE(2,80)'* OWENS & GIBBS (1964) REAERATION FORMULA IS USED '
      ELSEIF(IWQKA(1).EQ.4)THEN
        WRITE(2,80)'* MODIFIED OWENS & GIBBS REAERATION IS USED       '
      ENDIF
      IF(IWQVLIM.EQ.0)THEN
        WRITE(2,80)'* MACROALGAE GROWTH IS NOT LIMITED BY VELOCITY    '
      ELSEIF(IWQVLIM.EQ.1)THEN
        WRITE(2,80)'* MACROALGAE VELOCITY LIMIT, MICHAELIS-MENTON EQN.'
      ELSEIF(IWQVLIM.EQ.2)THEN
        WRITE(2,80)'*MACROALGAE VEL. LIMIT, 5-PARAM LOGISTIC FUNCTION'
      ENDIF
      WRITE(2,83)'* # OF ZONES FOR SPATIALLY VARYING PARAMETERS =',IWQZ
      IF(IWQZ.GT.NWQZ) STOP 'ERROR!! IWQZ SHOULD BE <= NWQZ'
      IF(IWQNC.EQ.1)THEN
        WRITE(2,80)'* WRITE NEGATIVE CONC. INFORMATION TO NEG-CONC.LOG'
      ELSE
        WRITE(2,80)'* NO WRTING OF NEGATIVE CONCENTRATION INFORMATION '
      ENDIF
      IF(IWQRST.EQ.1)THEN
        WRITE(2,80)'* WRITE SPATIAL DISTRIBUTIONS TO IWQORST          '
      ELSE
        WRITE(2,80)'* NO WRITING TO IWQORST                           '
      ENDIF
      WRITE(2,999)
      IF(IWQICI.EQ.1)THEN
        WRITE(2,80)'* SPATIALLY/TEMPORALLY-VARYING ICS FROM INWQICI   '
      ELSEIF(IWQICI.EQ.2)THEN
        WRITE(2,80)'* SPATIALLY/TEMPORALLY-VARYING ICS FROM INWQRST   '
      ELSEIF(IWQICI.EQ.3)THEN
        WRITE(2,80)'* SPATIALLY/TEMPORALLY CONSTANT FROM FIRST TIME SERIES   '
      ELSE
        WRITE(2,80)'* SPATIALLY/TEMPORALLY CONSTANT C44 INITIAL CONDITIONS'
      ENDIF
      IF(IWQAGR.EQ.1)THEN
        WRITE(2,80)'* SPATIALLY A/O TEMPORALLY-VARYING ALGAL KINETICS '
      ELSE
        WRITE(2,80)'* SPATIALLY/TEMPORALLY CONSTANT ALGAL KINETICS    '
      ENDIF
      IF(IWQSTL.EQ.1)THEN
        WRITE(2,80)'* SPATIALLY AND/OR TEMPORALLY-VARYING SETTLING VEL'
      ELSE
        WRITE(2,80)'* SPATIALLY/TEMPORALLY CONSTANT SETTLING VELOCITY '
      ENDIF
      IF(IWQSUN.GE.1)THEN
        WRITE(2,80)'* TEMPORALLY-VARYING IO & FD                      '
      ELSE
        WRITE(2,80)'* TEMPORALLY CONSTANT IO & FD                     '
      ENDIF
      IF(IWQNPL.EQ.1)THEN
        WRITE(2,80)'* SPATIALLY AND/OR TEMPORALLY-VARYING NPS INPUT   '
      ELSE
        WRITE(2,80)'* SPATIALLY/TEMPORALLY CONSTANT NPS INPUT         '
      ENDIF
      IF(IWQKIN.EQ.1)THEN
        WRITE(2,80)'* SPATIALLY VARYING KINETICS FROM KINETICS.INP    '
      ELSE
        WRITE(2,80)'* FILE KINETICS.INP NOT USED                      '
      ENDIF
      WRITE(2,999)
C *** C06
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      READ(1,*) IWQTS,TWQTSB,TWQTSE,WQTSDT, ISWQAVG, ISWQMIN, ISWQMAX,ISCOMP
      WRITE(2,*) IWQTS,TWQTSB,TWQTSE,WQTSDT, ISWQAVG, ISWQMIN, ISWQMAX,ISCOMP
C
C ISWQAVG > 0 TURNS ON BINARY FILE OUTPUT FOR WQ DAILY AVERAGES
C ISWQMIN > 0 TURNS ON BINARY FILE OUTPUT FOR WQ DAILY MINIMUMS
C ISWQMAX > 0 TURNS ON BINARY FILE OUTPUT FOR WQ DAILY MAXIMUMS
C ISCOMP  > 0 TURNS ON BINARY FILE OUTPUT FOR DO COMPONENT ANALYSIS
C
      IF(ISCOMP .GT. 0)THEN
        CALL WQZERO3
        CALL INITBIN3
      ENDIF
      IF(IWQTS.GT.NWQTS)THEN
        WRITE(2,80)'** IWQTS SHOULD BE <= NWQTS **                    '
        IWQTS=NWQTS
      ENDIF
      WRITE(2,84)
     &    '* TIME-SERIES OUTPUT FROM ', TWQTSB, ' DAY ',
     &    '                       TO ', TWQTSE, ' DAY ',
     &    '                    EVERY ', WQTSDT, ' HOUR',
     &    '                       AT ', IWQTS,  ' LOCATIONS',
     &    ' BIN FILE SWITCH ISWQAVG =', ISWQAVG,' (0=OFF)  ',
     &    ' BIN FILE SWITCH ISWQMIN =', ISWQMIN,' (0=OFF)  ',
     &    ' BIN FILE SWITCH ISWQMAX =', ISWQMAX,' (0=OFF)  ',
     &    ' BIN FILE SWITCH ISCOMP  =', ISCOMP, ' (0=OFF)  '
      WRITE(2,999)
C *** C07
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      DO M=1,2
        READ(1,90) TITLE(M)
      ENDDO
      IF(IWQTS.GE.1)THEN
        WRITE(2,80)': ICWQTS(I)=1, TIME-SERIES OUTPUT FOR VARIABLE I  '
        WRITE(2,80)': ICWQTS(I)\=1, NO TIME-SERIES OUTPUT FOR VAR. I  '
        WRITE(2,999)
        DO M=1,2
          WRITE(2,90) TITLE(M)
        ENDDO
        MM=0
        DO M=1,IWQTS
          READ(1,*) II,JJ,(ICWQTS(NW,M),NW=1,13)
          WRITE(2,*) II,JJ,(ICWQTS(NW,M),NW=1,13)
          READ(1,*) (ICWQTS(NW,M),NW=14,NTSWQV),ICWQTS(IDNOTRVA,M)
          WRITE(2,*) (ICWQTS(NW,M),NW=14,NTSWQV),ICWQTS(IDNOTRVA,M)
          I = XLOC(II)
          J = YLOC(JJ)
          L = LIJ(I,J)
          IF (CELL_INSIDE_DOMAIN(L)) THEN
            IF(IJCT(I,J).LT.1 .OR. IJCT(I,J).GT.8)THEN
              WRITE(2,86)  II,JJ,I,J,M
              WRITE(2,80)'ERROR!! INVALID (I,J): TIME-SERIES LOCATION'
              STOP
            END IF
            MM=MM+1
            LWQTS(MM)=LIJ(I,J)
            WRITE(2,94) II,JJ,(ICWQTS(NW,M),NW=1,NTSWQV+1)
          ENDIF
        ENDDO
      ENDIF
      IWQTSB = NINT(TWQTSB/DTD)
      IWQTSE = NINT(TWQTSE/DTD)
      IWQTSDT = NINT(WQTSDT*3600.0/DT)
      WRITE(2,999)
      WRITE(2,83)': TIME-SERIES STARTING TIME STEP (IN DT UNIT) =',IWQTSB
      WRITE(2,83)': TIME-SERIES ENDING TIME STEP (IN DT UNIT)   =',IWQTSE
      WRITE(2,83)': FREQUENCY OF TS OUTPUT  (IN DT UNIT)        =',IWQTSDT
C PMC      IF(MOD(IWQTSDT,IWQDT).NE.0)
C PMC     &    STOP 'ERROR!! IWQTSDT SHOULD BE MULTIPLE OF IWQDT'
  999 FORMAT(1X)
   90 FORMAT(A79)
   91 FORMAT(10I8)
   92 FORMAT(10F8.4)
   93 FORMAT(I8,3F8.4)
   94 FORMAT(2I5, 13I5, /, 10X, 9I5)
   95 FORMAT(A254)
   80 FORMAT(A50)
   81 FORMAT(A27, 4(F8.4,2X))
   82 FORMAT((A45, F8.3))
   83 FORMAT(A47, I10)
   84 FORMAT(3(A26,F10.4,A5,/), 5(A26,I8,A10,/))
   86 FORMAT(' I,J,M = ',3I10)
C
C CONSTANT PARAMETERS FOR ALGAE (SEE TABLE 3-1)
C8
C
      WRITE(2,999)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      READ(1,90) TITLE(1)
      WRITE(2,90) TITLE(1)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)

	IF(ISTRWQ(22) .GT. 0)THEN
        READ(1,*) WQKHNC,WQKHND,WQKHNG,WQKHNM,WQKHPC,WQKHPD,WQKHPG,		!AA Added read of CO2 half-saturation consts
     &    WQKHPM,WQKHS,WQSTOX,WQKHCO2C,WQKHCO2D,WQKHCO2G,WQKHCO2M
        WRITE(2,*) WQKHNC,WQKHND,WQKHNG,WQKHNM,WQKHPC,WQKHPD,WQKHPG,
     &    WQKHPM,WQKHS,WQSTOX,WQKHCO2C,WQKHCO2D,WQKHCO2G,WQKHCO2M
	  WRITE(2,80)'* HALF-SAT. CONSTANT (G/M^3) FOR NUTRIENT UPTAKE '
        WRITE(2,81)' : (KHNC, KHPC, KHCO2C)       = ', WQKHNC,WQKHPC,WQKHCO2C
        WRITE(2,81)' : (KHND, KHPD, KHS, KHCO2D)  = ', WQKHND,WQKHPD,WQKHS,WQKHCO2D
        WRITE(2,81)' : (KHNG, KHPG, KHCO2G)       = ', WQKHND,WQKHPG,WQKHCO2G
        WRITE(2,81)' : (KHNM, KHPM, KHCO2M)       = ', WQKHNM,WQKHPM,WQKHCO2M
      ELSE
	  READ(1,*) WQKHNC,WQKHND,WQKHNG,WQKHNM,WQKHPC,WQKHPD,WQKHPG,		!AA True if not using CO2 limitation
     &    WQKHPM,WQKHS,WQSTOX
        WRITE(2,*) WQKHNC,WQKHND,WQKHNG,WQKHNM,WQKHPC,WQKHPD,WQKHPG,
     &    WQKHPM,WQKHS,WQSTOX
	  WRITE(2,80)'* HALF-SAT. CONSTANT (G/M^3) FOR NUTRIENT UPTAKE '
        WRITE(2,81)' : (KHNC, KHPC)       = ', WQKHNC,WQKHPC
        WRITE(2,81)' : (KHND, KHPD, KHS)  = ', WQKHND,WQKHPD,WQKHS
        WRITE(2,81)' : (KHNG, KHPG)       = ', WQKHND,WQKHPG
        WRITE(2,81)' : (KHNM, KHPM)       = ', WQKHNM,WQKHPM
      ENDIF

      WRITE(2,82)'* SAL. WHERE MICROSYSTIS GROWTH IS HALVED  = ',WQSTOX
      WQSTOX = WQSTOX*WQSTOX
C
C9
C
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      ! *** BEGIN DSLLC BLOCK
      READ(1,95)LINE
      READ(LINE,*,END=109,ERR=109) WQKETSS,WQKECHL,WQCHLC,WQCHLD,WQCHLG,
     &    WQCHLM,WQDOPC,WQDOPD,WQDOPG,WQDOPM(1),WQKEPOM
  109 WRITE(2,*) WQKETSS,WQKECHL,WQCHLC,WQCHLD,WQCHLG,WQCHLM,WQDOPC,
     &    WQDOPD,WQDOPG, WQDOPM(1),WQKEPOM
      IF(ISTRAN(6).EQ.0)THEN
        WQKETSS=0.0
        WRITE(2,80)': SINCE TSS IS NOT MODELED, KETSS IS FORCED TO 0  '
      ENDIF
      WRITE(2,80)'* LIGHT EXTINC. COEFF. DUE TO TSS, CHL & POM      '
      WRITE(2,81)' : KETSS (/M PER G/M^3)  = ', WQKETSS
      WRITE(2,81)' : KECHL (/M PER MG/M^3) = ', WQKECHL
      IF(WQKECHL .LT. 0.0)THEN
        WRITE(2,80) '* USE RILEY (1956) EQUATION FOR WQKECHL          '
        WRITE(2,80) ' : KECHL = 0.054*CHL**0.667 + 0.0088*CHL         '
      ENDIF
      WRITE(2,81)' : KEPOM (/M PER G/M^3)  = ', WQKEPOM
      ! *** END DSLLC BLOCK
      WRITE(2,80)'* CARBON-TO-CHL RATIO (G C PER MG CHL)            '
      WRITE(2,81)' : (CCHLC, CCHLD, CCHLG) = ', WQCHLC,WQCHLD,WQCHLG
      WRITE(2,80)'* DEPTH (M) OF MAXIMUM ALGAL GROWTH               '
      WRITE(2,81)' : (DOPTC, DOPTD, DOPTG) = ', WQDOPC,WQDOPD,WQDOPG
      WQCHLC=1.0/(WQCHLC + 1.E-12)
      WQCHLD=1.0/(WQCHLD + 1.E-12)
      WQCHLG=1.0/(WQCHLG + 1.E-12)
      WQCHLM=1.0/(WQCHLM + 1.E-12)
C *** C10
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      READ(1,*) WQI0,WQISMIN,WQFD,WQCIA,WQCIB,WQCIC,WQCIM,REAC(1),
     &    PARADJ
      WRITE(2,*) WQI0,WQISMIN,WQFD,WQCIA,WQCIB,WQCIC,WQCIM,REAC(1),
     &    PARADJ
      WRITE(2,82)'*INITIAL IO (LY/D) AT WATER SURFACE       = ',WQI0
     &    ,' MINIMUM OPTIMUM SOLAR RADIATION (LY/D)   = ',WQISMIN
     &    ,' FRACTIONAL DAYLENGTH                     = ',WQFD
     &    ,' WEIGHTING FACTOR FOR RAD. AT CURRENT DAY = ',WQCIA
     &    ,' WEIGHTING FACTOR FOR RAD. AT (-1) DAY    = ',WQCIB
     &    ,' WEIGHTING FACTOR FOR RAD. AT (-2) DAYS   = ',WQCIC
     &    ,' FRACTION OF SOLAR RADIATION THAT IS PAR  = ',PARADJ
      WQI0=PARADJ*WQI0  !/(WQFD+1.E-18)  ! *** APPLY CONVERSION TO OPTIMAL LIGHT
      WQI1 = WQI0
      WQI2 = WQI0
      WQI3 = WQI0
      WQI0OPT = 0.0
      IF(IWQSUN .EQ. 2)WQISMIN = 0.0
C *** C11
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      READ(1,*)WQTMC1,WQTMC2,WQTMD1,WQTMD2,WQTMG1,WQTMG2,WQTMM1,WQTMM2,
     &    WQTMP1, WQTMP2
      WRITE(2,*)WQTMC1,WQTMC2,WQTMD1,WQTMD2,WQTMG1,WQTMG2,WQTMM1,
     &    WQTMM2,WQTMP1, WQTMP2
C *** C12
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      READ(1,*)WQKG1C,WQKG2C,WQKG1D,WQKG2D,WQKG1G,WQKG2G,WQKG1M,WQKG2M,
     &    WQKG1P, WQKG2P
      WRITE(2,*)WQKG1C,WQKG2C,WQKG1D,WQKG2D,WQKG1G,WQKG2G,WQKG1M,
     &    WQKG2M,WQKG1P, WQKG2P
      WRITE(2,80)'* LOWER OPTIMUM TEMP FOR ALGAL GROWTH (DEGC)     '
      WRITE(2,81)' : (TMC1, TMD1, TMG1   ) = ', WQTMC1,WQTMD1,WQTMG1
      WRITE(2,80)'* UPPER OPTIMUM TEMP FOR ALGAL GROWTH (DEGC)     '
      WRITE(2,81)' : (TMC2, TMD2, TMG2   ) = ', WQTMC2,WQTMD2,WQTMG2
C
C *** C13 CONSTANT PARAMETERS FOR ALGAE (SEE TABLE 3-1)
C
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      READ(1,*) WQTRC,WQTRD,WQTRG,WQTRM,WQKTBC,WQKTBD,WQKTBG,WQKTBM
      WRITE(2,*) WQTRC,WQTRD,WQTRG,WQTRM,WQKTBC,WQKTBD,WQKTBG,WQKTBM
      WRITE(2,80)'* REFERENCE TEMPERATURE FOR ALGAL METABOLISM (OC) '
      WRITE(2,81)' : (TRC, TRD, TRG)       = ', WQTRC,WQTRD,WQTRG
      WRITE(2,80)'* TEMPERATURE EFFECT FOR ALGAL METABOLISM         '
      WRITE(2,81)' : (KTBC, KTBD, KTBG)    = ', WQKTBC,WQKTBD,WQKTBG
C
C *** C14 CONSTANT PARAMETERS FOR CARBON (SEE TABLE 3-2)
C
      WRITE(2,999)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      READ(1,90) TITLE(1)
      WRITE(2,90) TITLE(1)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      READ(1,*) WQFCRP,WQFCLP,WQFCDP,WQFCDC,WQFCDD,WQFCDG,WQKHRC,WQKHRD,WQKHRG
      WRITE(2,*) WQFCRP,WQFCLP,WQFCDP,WQFCDC,WQFCDD,WQFCDG,WQKHRC,WQKHRD,WQKHRG
      WRITE(2,80)'* CARBON DISTRIBUTION COEFF FOR ALGAL PREDATION   '
      WRITE(2,81)' : (FCRP, FCLP, FCDP)    = ', WQFCRP,WQFCLP,WQFCDP
      WRITE(2,80)'* CARBON DISTRIBUTION COEFF FOR ALGAL METABOLISM  '
      WRITE(2,81)' : (FCDC, FCDD, FCDG)    = ', WQFCDC,WQFCDD,WQFCDG
      WRITE(2,80)'* HALF-SAT. CONSTANT (GO/M*3) FOR ALGAL DOC EXCRET'
      WRITE(2,81)' : (KHRC, KHRD, KHRG)    = ', WQKHRC,WQKHRD,WQKHRG
      CFCDCWQ = 1.0 - WQFCDC
      CFCDDWQ = 1.0 - WQFCDD
      CFCDGWQ = 1.0 - WQFCDG
      XC = ABS(1.0 - (WQFCRP+WQFCLP+WQFCDP))
      IF(XC .GT. 0.0001)THEN
        WRITE(2,*)
        WRITE(2,*) ' WARNING!  FCRP+FCLP+FCDP NOT EQUAL TO 1.0'
        WRITE(2,*)
      ENDIF
C *** C15
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      READ(1,90) TITLE(1)
      WRITE(2,90) TITLE(1)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      READ(1,*) WQFCRPM,WQFCLPM,WQFCDPM,WQFCDM, WQKHRM(1)
      WRITE(2,*) WQFCRPM,WQFCLPM,WQFCDPM,WQFCDM, WQKHRM(1)
C *** C16
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      READ(1,*)WQKRC,WQKLC,WQKDC(1),WQKRCALG,WQKLCALG,WQKDCALG,WQKDCALM(1)
      WRITE(2,*)WQKRC,WQKLC,WQKDC(1),WQKRCALG,WQKLCALG,WQKDCALG,WQKDCALM(1)
      WRITE(2,80)'* MINIMUM DISSOLUTION RATE (/DAY) OF ORGANIC C    '
      WRITE(2,81)' : (KRC, KLC, KDC)       = ', WQKRC,WQKLC,WQKDC(1)
      WRITE(2,80)'* CONSTANT RELATING DISSOLUTION RATE TO ALGAE     '
      WRITE(2,81)' : (KRCALG,KLCALG,KDCALG)= ', WQKRCALG,WQKLCALG,WQKDCALG
C *** C17
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      READ(1,*) WQTRHDR,WQTRMNL,WQKTHDR,WQKTMNL,WQKHORDO,WQKHDNN,WQAANOX
      WRITE(2,*) WQTRHDR,WQTRMNL,WQKTHDR,WQKTMNL,WQKHORDO,WQKHDNN,WQAANOX
      WRITE(2,80)'* REFERENCE TEMP FOR HYDROLYSIS/MINERALIZATION(OC)'
      WRITE(2,81)' : (TRHDR, TRMNL)        = ', WQTRHDR,WQTRMNL
      WRITE(2,80)'* TEMPERATURE EFFECT ON HYDROLYSIS/MINERALIZATION '
      WRITE(2,81)' : (KTHDR, KTMNL)        = ', WQKTHDR,WQKTMNL
      WRITE(2,80)'* HALF-SAT. CONSTANT FOR OXIC RESP/DENITRIFICATION'
      WRITE(2,81)' : (KHORDO, KHDNN)       = ', WQKHORDO,WQKHDNN
      WRITE(2,80)'* RATION OF DENITRIFICATION TO OXIC DOC RESP      '
      WRITE(2,81)' : (AANOX)               = ', WQAANOX
      WQAANOX = WQAANOX*WQKHORDO
C
C *** C18 CONSTANT PARAMETERS FOR PHOSPHORUS (TABLE 3-3)
C
      WRITE(2,999)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      READ(1,90) TITLE(1)
      WRITE(2,90) TITLE(1)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      READ(1,*) WQFPRP,WQFPLP,WQFPDP,WQFPIP,WQFPRC,WQFPRD,WQFPRG,
     &    WQFPLC,WQFPLD,WQFPLG
      WRITE(2,*) WQFPRP,WQFPLP,WQFPDP,WQFPIP,WQFPRC,WQFPRD,WQFPRG,
     &    WQFPLC,WQFPLD,WQFPLG
      WRITE(2,80)'* PHOSPHORUS DISTRIBUTION COEF FOR ALGAL PREDATION'
      WRITE(2,81)' : (FPRP,FPLP,FPDP,FPIP) = ', WQFPRP,WQFPLP,WQFPDP,WQFPIP
      WRITE(2,80)'* PHOSPHORUS DIST COEF OF RPOP FOR ALGAL METABOLIS'
      WRITE(2,81)' : (FPRC, FPRD, FPRG)    = ', WQFPRC,WQFPRD,WQFPRG
      WRITE(2,80)'* PHOSPHORUS DIST COEF OF LPOP FOR ALGAL METABOLIS'
      WRITE(2,81)' : (FPLC, FPLD, FPLG)    = ', WQFPLC,WQFPLD,WQFPLG
      XP = ABS(1.0 - (WQFPRP+WQFPLP+WQFPDP+WQFPIP))
      IF(XP .GT. 0.0001)THEN
        WRITE(2,*)
        WRITE(2,*) ' WARNING!  FPRP+FPLP+FPDP+FPIP NOT EQUAL TO 1.0'
        WRITE(2,*)
      ENDIF
C *** C19
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      READ(1,90) TITLE(1)
      WRITE(2,90) TITLE(1)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      IF(ISSKIP .EQ. 0) READ(1,999)
      READ(1,*)WQFPRPM,WQFPLPM,WQFPDPM,WQFPIPM,WQFPRM,WQFPLM,WQAPCM
      WRITE(2,*)WQFPRPM,WQFPLPM,WQFPDPM,WQFPIPM,WQFPRM,WQFPLM,WQAPCM
C *** C20
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      READ(1,*) WQFPDC,WQFPDD,WQFPDG,WQFPDM,WQFPIC,WQFPID,WQFPIG,WQFPIM,WQKPO4P
      WRITE(2,*) WQFPDC,WQFPDD,WQFPDG,WQFPDM,WQFPIC,WQFPID,WQFPIG,WQFPIM,WQKPO4P
      IF(IWQSRP.NE.1 .AND. IWQSRP.NE.2)THEN
        WQKPO4P = 0.0
        WRITE(2,80)': NO SORPTION OF PO4T/SA, SO KPO4P IS FORCED TO 0 '
      ENDIF
      WRITE(2,80)'* PHOSPHORUS DIST COEF OF DOP FOR ALGAL METABOLISM'
      WRITE(2,81)' : (FPDC, FPDD, FPDG)    = ', WQFPDC,WQFPDD,WQFPDG
      WRITE(2,80)'* PHOSPHORUS DIST COEF OF NH4 FOR ALGAL METABOLISM'
      WRITE(2,81)' : (FPIC, FPID, FPIG)    = ', WQFPIC,WQFPID,WQFPIG
      WRITE(2,82)'* PARITITION COEFF FOR SORBED/DISSOLVED PO4 =',WQKPO4P
      XPC = ABS(1.0 - (WQFPRC+WQFPLC+WQFPDC+WQFPIC))
      IF(XPC .GT. 0.0001)THEN
        WRITE(2,*)
        WRITE(2,*) ' WARNING!  FPRC+FPLC+FPDC+FPIC NOT EQUAL TO 1.0'
        WRITE(2,*)
      ENDIF
      XPD = ABS(1.0 - (WQFPRD+WQFPLD+WQFPDD+WQFPID))
      IF(XPD .GT. 0.0001)THEN
        WRITE(2,*)
        WRITE(2,*) ' WARNING!  FPRD+FPLD+FPDD+FPID NOT EQUAL TO 1.0'
        WRITE(2,*)
      ENDIF
      XPG = ABS(1.0 - (WQFPRG+WQFPLG+WQFPDG+WQFPIG))
      IF(XPG .GT. 0.0001)THEN
        WRITE(2,*)
        WRITE(2,*) ' WARNING!  FPRG+FPLG+FPDG+FPIG NOT EQUAL TO 1.0'
        WRITE(2,*)
      ENDIF
C *** C21
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      READ(1,*) WQKRP,WQKLP,WQKDP,WQKRPALG,WQKLPALG,WQKDPALG,WQCP1PRM,WQCP2PRM,WQCP3PRM
      WRITE(2,*) WQKRP,WQKLP,WQKDP,WQKRPALG,WQKLPALG,WQKDPALG,WQCP1PRM,WQCP2PRM,WQCP3PRM
      WRITE(2,80)'* MINIMUM HYDROLYSIS RATE (/DAY) OF ORGANIC P     '
      WRITE(2,81)' : (KRP, KLP, KDP)       = ', WQKRP,WQKLP,WQKDP
      WRITE(2,80)'* CONSTANT RELATING HYDROLYSIS RATE TO ALGAE      '
      WRITE(2,81)' : (KRPALG,KLPALG,KDPALG)= ', WQKRPALG,WQKLPALG,WQKDPALG
      WRITE(2,80)'* CONSTANT USED IN DETERMINING P-TO-C RATIO       '
      WRITE(2,81)' : (CPPRM1,CPPRM2,CPPRM3)= ', WQCP1PRM,WQCP2PRM,WQCP3PRM
C
C *** C22 CONSTANT PARAMETERS FOR NITROGEN (TABLE 3-4)
C
      WRITE(2,999)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      READ(1,90) TITLE(1)
      WRITE(2,90) TITLE(1)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      READ(1,*) WQFNRP,WQFNLP,WQFNDP,WQFNIP,WQFNRC,WQFNRD,WQFNRG,
     &    WQFNLC,WQFNLD,WQFNLG
      WRITE(2,*) WQFNRP,WQFNLP,WQFNDP,WQFNIP,WQFNRC,WQFNRD,WQFNRG,
     &    WQFNLC,WQFNLD,WQFNLG
      WRITE(2,80)'* NITROGEN DISTRIBUTION COEFF FOR ALGAL PREDATION '
      WRITE(2,81)' : (FNRP,FNLP,FNDP,FNIP) = ', WQFNRP,WQFNLP,WQFNDP,
     &    WQFNIP
      WRITE(2,80)'* NITROGEN DIST COEF OF RPON FOR ALGAL METABOLISM '
      WRITE(2,81)' : (FNRC, FNRD, FNRG)    = ', WQFNRC,WQFNRD,WQFNRG
      WRITE(2,80)'* NITROGEN DIST COEF OF LPON FOR ALGAL METABOLISM '
      WRITE(2,81)' : (FNLC, FNLD, FNLG)    = ', WQFNLC,WQFNLD,WQFNLG
      XN = ABS(1.0 - (WQFNRP+WQFNLP+WQFNDP+WQFNIP))
      IF(XN .GT. 0.0001)THEN
        WRITE(2,*)
        WRITE(2,*) ' WARNING!  FNRP+FNLP+FNDP+FNIP NOT EQUAL TO 1.0'
        WRITE(2,*)
      ENDIF
      WRITE(2,999)
C *** C23
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      READ(1,90) TITLE(1)
      WRITE(2,90) TITLE(1)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      READ(1,*)WQFNRPM,WQFNLPM,WQFNDPM,WQFNIPM,WQFNRM,WQFNLM
      WRITE(2,*)WQFNRPM,WQFNLPM,WQFNDPM,WQFNIPM,WQFNRM,WQFNLM
C *** C24
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      READ(1,*) WQFNDC,WQFNDD,WQFNDG,WQFNDM,WQFNIC,WQFNID,WQFNIG,
     &    WQFNIM,WQANCC,WQANCD,WQANCG,WQANCM
      WRITE(2,*) WQFNDC,WQFNDD,WQFNDG,WQFNDM,WQFNIC,WQFNID,WQFNIG,
     &    WQFNIM,WQANCC,WQANCD,WQANCG,WQANCM
      WRITE(2,80)'* NITROGEN DIST COEF OF DON FOR ALGAL METABOLISM  '
      WRITE(2,81)' : (FNDC, FNDD, FNDG)    = ', WQFNDC,WQFNDD,WQFNDG
      WRITE(2,80)'* NITROGEN DIST COEF OF NH4 FOR ALGAL METABOLISM  '
      WRITE(2,81)' : (FNIC, FNID, FNIG)    = ', WQFNIC,WQFNID,WQFNIG
      WRITE(2,80)'* NITROGEN-TO-CARBON RATIO IN ALGAE               '
      WRITE(2,81)' : (ANCC, ANCD, ANCG)    = ', WQANCC,WQANCD,WQANCG
      XNC = ABS(1.0 - (WQFNRC+WQFNLC+WQFNDC+WQFNIC))
      IF(XNC .GT. 0.0001)THEN
        WRITE(2,*)
        WRITE(2,*) ' WARNING!  FNRC+FNLC+FNDC+FNIC NOT EQUAL TO 1.0'
        WRITE(2,*)
      ENDIF
      XND = ABS(1.0 - (WQFNRD+WQFNLD+WQFNDD+WQFNID))
      IF(XND .GT. 0.0001)THEN
        WRITE(2,*)
        WRITE(2,*) ' WARNING!  FNRD+FNLD+FNDD+FNID NOT EQUAL TO 1.0'
        WRITE(2,*)
      ENDIF
      XNG = ABS(1.0 - (WQFNRG+WQFNLG+WQFNDG+WQFNIG))
      IF(XNG .GT. 0.0001)THEN
        WRITE(2,*)
        WRITE(2,*) ' WARNING!  FNRG+FNLG+FNDG+FNIG NOT EQUAL TO 1.0'
        WRITE(2,*)
      ENDIF
C *** C25
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      READ(1,*) WQANDC,WQNITM,WQKHNDO,WQKHNN,WQTNIT,WQKN1,WQKN2
      WRITE(2,*) WQANDC,WQNITM,WQKHNDO,WQKHNN,WQTNIT,WQKN1,WQKN2
      WRITE(2,82)'* MASS NO3 REDUCED PER DOC OXIDIZED (GN/GC)= ',WQANDC
     &    ,'* MAXIMUM NITRIFICATION RATE (G N /M^3/D)  = ',WQNITM
     &    ,'  REFERENCE TEMP FOR NITRIFICATION (DEGC)  = ',WQTNIT
      WRITE(2,80)'* NITRIFICATION HALF-SAT CONSTANT FOR DO & NH4    '
      WRITE(2,81)' : (KHNITDO, KHNITN)     = ', WQKHNDO,WQKHNN
      WRITE(2,80)'* SUB & SUPER-OPTIMUM TEMP EFFECT ON NITRIFICATION'
      WRITE(2,81)' : (KNIT1, KNIT2)        = ', WQKN1,WQKN2
C *** C26
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      READ(1,*) WQKRN,WQKLN,WQKDN,WQKRNALG,WQKLNALG,WQKDNALG
      WRITE(2,*) WQKRN,WQKLN,WQKDN,WQKRNALG,WQKLNALG,WQKDNALG
      WRITE(2,80)'* MINIMUM HYDROLYSIS RATE (/DAY) OF ORGANIC N     '
      WRITE(2,81)' : (KRN, KLN, KDN)       = ', WQKRN,WQKLN,WQKDN
      WRITE(2,80)'* CONSTANT RELATING HYDROLYSIS RATE TO ALGAE      '
      WRITE(2,81)' : (KRNALG,KLNALG,KDNALG)= ', WQKRNALG,WQKLNALG,WQKDNALG
C
C *** C27 CONSTANT PARAMETERS FOR SILICA (TABLE 3-5)
C
      WRITE(2,999)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      READ(1,90) TITLE(1)
      WRITE(2,90) TITLE(1)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      READ(1,*) WQFSPP,WQFSIP,WQFSPD,WQFSID,WQASCD,WQKSAP,WQKSU,
     &    WQTRSUA,WQKTSUA
      WRITE(2,*) WQFSPP,WQFSIP,WQFSPD,WQFSID,WQASCD,WQKSAP,WQKSU,
     &    WQTRSUA,WQKTSUA
      IF(IWQSRP.NE.1 .AND. IWQSRP.NE.2)THEN
        WQKSAP = 0.0
        WRITE(2,80)': NO SORPTION OF PO4T/SA, SO KSAP IS FORCED TO 0  '
      ENDIF
      WRITE(2,80)'* SILICA DISTRIBUTION COEFF FOR DIATOM PREDATION  '
      WRITE(2,81)' : (FSPP, FSIP)          = ', WQFSPP,WQFSIP
      WRITE(2,80)'* SILICA DISTRIBUTION COEFF FOR DIATOM METABOLISM '
      WRITE(2,81)' : (FSPD, FSID)          = ', WQFSPD,WQFSID
      WRITE(2,82)'*SILICA-TO-CARBON RATIO IN DIATOMS        = ',WQASCD
     &    ,'*PARITITION COEFF FOR SORBED/DISSOLVED SA = ',WQKSAP
     &    ,'*DISSOLUTION RATE (/D) OF PSI             = ',WQKSU
     &    ,' REFERENCE TEMP FOR PSI DISSOLUTION (OC)  = ',WQTRSUA
     &    ,' TEMPERATURE EFFECT ON PSI DISSOLUTION    = ',WQKTSUA
C
C *** C28 CONSTANT PARAMETERS FOR COD & DO (TABLE 3-6)
C
      WRITE(2,999)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      READ(1,90) TITLE(1)
      WRITE(2,90) TITLE(1)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      READ(1,*) WQAOCR,WQAONT, WQKRO(1), WQKTR(1),WQKHCOD(1),WQKCD(1),
     &    WQTRCOD, WQKTCOD, WQAOCRPM, WQAOCRRM
      WRITE(2,*) WQAOCR,WQAONT, WQKRO(1), WQKTR(1),WQKHCOD(1),WQKCD(1),
     &    WQTRCOD, WQKTCOD, WQAOCRPM, WQAOCRRM
      WRITE(2,82)'* DO-TO-CARBON RATIO IN RESPIRATION        = ',WQAOCR
     &    ,':MASS DO CONSUMED PER NH4-N NITRIFIED     = ',WQAONT
     &    ,':PROPORN. CONSTANT FOR DO-REAERATION (MKS)= ',WQKRO(1)
     &    ,' TEMPERATURE EFFECT ON DO-REAERATION      = ',WQKTR(1)
     &    ,'*HALF-SAT CONSTANT OF DO FOR COD (GO2/M^3)= ',WQKHCOD(1)
     &    ,':OXIDATION RATE OF COD (/DAY)             = ',WQKCD(1)
     &    ,'  REFERENCE TEMP FOR COD OXIDATION (OC)    = ',WQTRCOD
     &    ,'  TEMPERATURE EFFECT ON COD OXIDATION      = ',WQKTCOD
     &    ,': DO-TO-CARBON RATIO MACROALGAE PHOTOSYNTH = ',WQAOCRPM
     &    ,': DO-TO-CARBON RATIO MACROALGAE RESPIRATION= ',WQAOCRRM
C
C *** C29 CONSTANT PARAMETERS FOR TAM & FCB (TABLE 3-7)
C
      WRITE(2,999)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      READ(1,90) TITLE(1)
      WRITE(2,90) TITLE(1)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      READ(1,*) WQKHBMF,WQBFTAM,WQTTAM,WQKTAM,WQTAMDMX,WQKDOTAM,
     &    WQKFCB,WQTFCB
      WRITE(2,*) WQKHBMF,WQBFTAM,WQTTAM,WQKTAM,WQTAMDMX,WQKDOTAM,
     &    WQKFCB,WQTFCB
      WRITE(2,82)
     &    '* DO WHERE TAM RELEASE IS HALF ANOXIC ONE  = ',WQKHBMF
     &    ,'  ANOXIC RELEASE OF TAM (MOL/M^2/D)        = ',WQBFTAM
     &    ,'  REFERENCE TEMP FOR TAM RELEASE (OC)      = ',WQTTAM
     &    ,'  TEMPERATURE EFFECT ON TAM RELEASE        = ',WQKTAM
     &    ,': TAM SOLUBILITY AT ANOXIC COND. (MOL/M^3) = ',WQTAMDMX
     &    ,'  CONSTANT RELATING TAM SOLUBILITY TO DO   = ',WQKDOTAM
     &    ,'* FIRST-ORDER DIE-OFF RATE AT 20OC (/D)    = ',WQKFCB
     &    ,'  TEMPERATURE EFFECT ON BACTERIA DIE-OFF   = ',WQTFCB
C
C SET UP LOOK-UP TABLE FOR TEMPERATURE DEPENDENCY OVER -15 OC TO 40 OC
C
      IF(NWQTD.GT.1000)THEN
        STOP 'ERROR!! NWQTD MUST BE <= 1000'
      ENDIF
      WQTDMIN=-10
      WQTDMAX=50
      WTEMP=WQTDMIN
      WQTDINC=(WQTDMAX-WQTDMIN)/NWQTD
      DO M=1,NWQTD
C        WTEMP =0.25*REAL(M)-30.25
        WQTDTEMP(M)=WTEMP
        IF(WTEMP.LT.WQTMC1)THEN
          WQTDGC(M) = EXP( -WQKG1C*(WTEMP-WQTMC1)*(WTEMP-WQTMC1) )
        ELSEIF(WTEMP.GT.WQTMC2)THEN
          WQTDGC(M) = EXP( -WQKG2C*(WTEMP-WQTMC2)*(WTEMP-WQTMC2) )
        ELSE
          WQTDGC(M)=1.0
        ENDIF
        IF(WTEMP.LT.WQTMD1)THEN
          WQTDGD(M) = EXP( -WQKG1D*(WTEMP-WQTMD1)*(WTEMP-WQTMD1) )
        ELSEIF(WTEMP.GT.WQTMD2)THEN
          WQTDGD(M) = EXP( -WQKG2D*(WTEMP-WQTMD2)*(WTEMP-WQTMD2) )
        ELSE
          WQTDGD(M)=1.0
        ENDIF
        IF(WTEMP.LT.WQTMG1)THEN
          WQTDGG(M) = EXP( -WQKG1G*(WTEMP-WQTMG1)*(WTEMP-WQTMG1) )
        ELSEIF(WTEMP.GT.WQTMG2)THEN
          WQTDGG(M) = EXP( -WQKG2G*(WTEMP-WQTMG2)*(WTEMP-WQTMG2) )
        ELSE
          WQTDGG(M)=1.0
        ENDIF
        IF(IDNOTRVA.GT.0)THEN
          IF(WTEMP.LT.WQTMM1)THEN
            WQTDGM(M) = EXP( -WQKG1M*(WTEMP-WQTMM1)*(WTEMP-WQTMM1) )
          ELSEIF(WTEMP.GT.WQTMM2)THEN
            WQTDGM(M) = EXP( -WQKG2M*(WTEMP-WQTMM2)*(WTEMP-WQTMM2) )
          ELSE
            WQTDGM(M)=1.0
          ENDIF
          WQTDRM(M) = EXP( WQKTBM*(WTEMP-WQTRM) )
        ENDIF
C
C  THE FOLLOWING WQTDGP VARIABLE IS A TEMPERATURE RELATED ADJUSTMENT
C  TO THE PREDATION AND/OR BASAL MATABOLISM RATE TO ALLOW DIATOMS
C  TO BLOOM IN WINTER (OR OTHER TIME OF YEAR).
C
        WQTDGP(M)=1.
        IF(WTEMP.LT.WQTMP1)WQTDGP(M) = EXP(-WQKG1P*(WTEMP-WQTMP1)*(WTEMP-WQTMP1) )
        IF(WTEMP.GT.WQTMD2)WQTDGP(M) = EXP(-WQKG2P*(WTEMP-WQTMP2)*(WTEMP-WQTMP2) )
  555 FORMAT(F7.2,4E12.4)
        WQTDRC(M) = EXP( WQKTBC*(WTEMP-WQTRC) )
        WQTDRD(M) = EXP( WQKTBD*(WTEMP-WQTRD) )
        WQTDRG(M) = EXP( WQKTBG*(WTEMP-WQTRG) )
        WQTDHDR(M) = EXP( WQKTHDR*(WTEMP-WQTRHDR) )
        WQTDMNL(M) = EXP( WQKTMNL*(WTEMP-WQTRMNL) )
        IF(WTEMP.GT.WQTNIT)THEN
          WQTDNIT(M)=WQNITM*EXP( WQKN1*(WTEMP-WQTNIT)*(WQTNIT-WTEMP) )
        ELSE
          WQTDNIT(M)=WQNITM*EXP( WQKN2*(WTEMP-WQTNIT)*(WQTNIT-WTEMP) )
        ENDIF
        WQKSUA(M) = WQKSU * EXP( WQKTSUA*(WTEMP-WQTRSUA) )
        WQKCOD(M,1) = WQKCD(1) * EXP( WQKTCOD*(WTEMP-WQTRCOD) )
        TT20 = WTEMP-20.0
        WQTDKR(M,1) = WQKTR(1)**TT20
        WRITE(2,2222)M,WQKTR(1),WQTDKR(M,1)
        WQTDTAM(M) = WQKHBMF * WQBFTAM * EXP( WQKTAM*(WTEMP-WQTTAM) )
        WQTT = WQKFCB * WQTFCB**TT20 * DTWQO2
        WQTD1FCB(M) = 1.0 - WQTT
        WQTD2FCB(M) = 1.0 / (1.0 + WQTT)
        WTEMP=WTEMP + WQTDINC
      ENDDO
C
C  *** C30
C PARAMETERS FOR WATER-QUALITY STATE-VARIABLE TIME SERIES
C
      WRITE(2,999)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      READ(1,90) TITLE(1)
      WRITE(2,90) TITLE(1)
      DO M=1,5
        READ(1,90) TITLE(M)
        WRITE(2,90) TITLE(M)
      ENDDO
      READ(1,*) (NWQCSR(NW),NW=1,NWQV)
      WRITE(2,*) (NWQCSR(NW),NW=1,NWQV)
      WRITE(2,970)(NWQCSR(NW),NW=1,NWQV)
      ! *** SAVE THE NUMBER OF WQ TIME SERIES
      DO NW=1,NWQV
        NT=4+NTOX+NSED+NSND+NW
        NCSER(NT)=NWQCSR(NW)
      ENDDO
C
C  *** C31
C PARAMETERS FOR OPEN BOUNDARY CONDITIONS
C
      WRITE(2,999)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      READ(1,90) TITLE(1)
      WRITE(2,90) TITLE(1)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      READ(1,*) NWQOBS_GL,NWQOBW_GL,NWQOBE_GL,NWQOBN_GL
      WRITE(2,*) NWQOBS_GL,NWQOBW_GL,NWQOBE_GL,NWQOBN_GL
      WRITE(2,23)'* # OF OPEN BDRY CELLS ON SOUTH BDRY       = ',NWQOBS_GL
      WRITE(2,23)'* # OF OPEN BDRY CELLS ON WEST BDRY        = ',NWQOBW_GL
      WRITE(2,23)'* # OF OPEN BDRY CELLS ON EAST BDRY        = ',NWQOBE_GL
      WRITE(2,23)'* # OF OPEN BDRY CELLS ON NORTH BDRY       = ',NWQOBN_GL
      IF(NWQOBS_GL.GT.NBBSM) STOP 'ERROR!! NWQOBS SHOULD <= NBBSM'
      IF(NWQOBW_GL.GT.NBBWM) STOP 'ERROR!! NWQOBW SHOULD <= NBBWM'
      IF(NWQOBE_GL.GT.NBBEM) STOP 'ERROR!! NWQOBE SHOULD <= NBBEM'
      IF(NWQOBN_GL.GT.NBBNM) STOP 'ERROR!! NWQOBN SHOULD <= NBBNM'
      WRITE(2,999)
      WRITE(2,80)'* CONSTANT OBC AT (ICBX(M),JCBX(M)) IF IWQOBX(M)=0'
      WRITE(2,80)': READ TIME-SERIES OBCS IWQOBX TIMES IF IWQOBX > 0'
C
C  *** C32
C SOUTH BOUNDARY
C
      WRITE(2,999)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      READ(1,90) TITLE(1)
      WRITE(2,90) TITLE(1)
      DO M=1,5
        READ(1,90) TITLE(M)
        WRITE(2,90) TITLE(M)
      ENDDO
      NWQOBS=0
      M=0
      IF(NWQOBS_GL.GT.0)THEN
        DO MM=1,NWQOBS_GL
          READ(1,*) I_TEMP,J_TEMP,(IWQOBS_GL(MM,NW),NW=1,NWQV)
          I_TEMP=XLOC(I_TEMP)
          J_TEMP=YLOC(J_TEMP)
          L = LIJ(I_TEMP, J_TEMP)
          IF (CELL_INSIDE_DOMAIN(L)) THEN
            NWQOBS=NWQOBS+1
            M=M+1
            IWQCBS(M) = I_TEMP
            JWQCBS(M) = J_TEMP
            IWQOBS(M,1:NWQV)=IWQOBS_GL(MM,1:NWQV)
            ! *** CONCENTRATION ASSIGNMENTS
            IF(IWQCBS(M).EQ.ICBS(M).AND.JWQCBS(M).EQ.JCBS(M))THEN
              DO NW=1,NWQV
                NT=4+NTOX+NSED+NSND+NW
                NCSERS(M,NT)=IWQOBS(M,NW)
              ENDDO
            ELSE
              STOP '  WQ: SOUTH OBC: MISMATCH BETWEEN NCBS & NWQOBS'
            ENDIF
            WRITE(2,*) IWQCBS(M),JWQCBS(M),(IWQOBS(M,NW),NW=1,NWQV)
            WRITE(2,969) IWQCBS(M),JWQCBS(M),(IWQOBS(M,NW),NW=1,NWQV)
          ENDIF
        ENDDO
      ENDIF
C
C  *** C33
C: CONSTANT BOTTOM AND SURFACE OBCS
C
      WRITE(2,999)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      DO M=1,5
        READ(1,90) TITLE(M)
        WRITE(2,90) TITLE(M)
      ENDDO
      NWQOBS=0
      M=0
      IF(NWQOBS_GL.GT.0)THEN
        DO MM=1,NWQOBS_GL
          READ(1,*) I_TEMP,J_TEMP,(WQOBCS_GL(MM,1,NW),NW=1,NWQV)
          I_TEMP=XLOC(I_TEMP)
          J_TEMP=YLOC(J_TEMP)
          L = LIJ(I_TEMP, J_TEMP)
          IF (CELL_INSIDE_DOMAIN(L)) THEN
            NWQOBS=NWQOBS+1
            M=M+1
            IWQCBS(M) = I_TEMP
            JWQCBS(M) = J_TEMP
            WQOBCS(M,1,1:NWQV)=WQOBCS_GL(MM,1,1:NWQV)
            ! *** CONCENTRATION ASSIGNMENTS, BOTTOM
            IF(IWQCBS(M).EQ.ICBS(M).AND.JWQCBS(M).EQ.JCBS(M))THEN
              DO NW=1,NWQV
                NT=4+NTOX+NSED+NSND+NW
                CBS(M,1,NT)=WQOBCS(M,1,NW)
              ENDDO
            ENDIF
            WRITE(2,*) IWQCBS(M),JWQCBS(M),(WQOBCS(M,1,NW),NW=1,NWQV)
            WRITE(2,97) IWQCBS(M),JWQCBS(M),(WQOBCS(M,1,NW),NW=1,NWQV)
          ENDIF
        ENDDO
      ENDIF
C  *** C34
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      WRITE(2,999)
      DO M=1,5
        READ(1,90) TITLE(M)
        WRITE(2,90) TITLE(M)
      ENDDO
      NWQOBS=0
      M=0
      IF(NWQOBS_GL.GT.0)THEN
        DO MM=1,NWQOBS_GL
          READ(1,*) I_TEMP,J_TEMP,(WQOBCS_GL(MM,2,NW),NW=1,NWQV)
          I_TEMP=XLOC(I_TEMP)
          J_TEMP=YLOC(J_TEMP)
          L=LIJ(I_TEMP,J_TEMP)
          IF (CELL_INSIDE_DOMAIN(L)) THEN
            NWQOBS=NWQOBS+1
            M=M+1
            IWQCBS(M) = I_TEMP
            JWQCBS(M) = J_TEMP
            WQOBCS(M,2,1:NWQV)=WQOBCS_GL(MM,2,1:NWQV)
          ! *** CONCENTRATION ASSIGNMENTS, TOP
            IF(IWQCBS(M).EQ.ICBS(M).AND.JWQCBS(M).EQ.JCBS(M))THEN
              DO NW=1,NWQV
                NT=4+NTOX+NSED+NSND+NW
                CBS(M,2,NT)=WQOBCS(M,2,NW)
              ENDDO
            ENDIF
            WRITE(2,*) IWQCBS(M),JWQCBS(M),(WQOBCS(M,2,NW),NW=1,NWQV)
            WRITE(2,97) IWQCBS(M),JWQCBS(M),(WQOBCS(M,2,NW),NW=1,NWQV)
          ENDIF
        ENDDO
      ENDIF
C
C  *** C35
C WEST BOUNDARY
C
      WRITE(2,999)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      READ(1,90) TITLE(1)
      WRITE(2,90) TITLE(1)
      DO M=1,5
        READ(1,90) TITLE(M)
        WRITE(2,90) TITLE(M)
      ENDDO
      NWQOBS=0
      M=0
      IF(NWQOBW_GL.GT.0)THEN
        DO MM=1,NWQOBW_GL
          READ(1,*) I_TEMP,J_TEMP,(IWQOBW_GL(MM,NW),NW=1,NWQV)
          I_TEMP=XLOC(I_TEMP)
          J_TEMP=YLOC(J_TEMP)
          L=LIJ(I_TEMP,J_TEMP)
          IF (CELL_INSIDE_DOMAIN(L)) THEN
            NWQOBW=NWQOBW+1
            M=M+1
            IWQCBW(M) = I_TEMP
            JWQCBW(M) = J_TEMP
            IWQOBW(M,1:NWQV)=IWQOBW_GL(MM,1:NWQV)
          ! *** CONCENTRATION ASSIGNMENTS
            IF(IWQCBW(M).EQ.ICBW(M).AND.JWQCBW(M).EQ.JCBW(M))THEN
              DO NW=1,NWQV
                NT=4+NTOX+NSED+NSND+NW
                NCSERW(M,NT)=IWQOBW(M,NW)
              ENDDO
            ELSE
              STOP '  WQ: WST OBC: MISMATCH BETWEEN NCBW & NWQOBW'
            ENDIF
            WRITE(2,*) IWQCBW(M),JWQCBW(M),(IWQOBW(M,NW),NW=1,NWQV)
            WRITE(2,969) IWQCBW(M),JWQCBW(M),(IWQOBW(M,NW),NW=1,NWQV)
          ENDIF
        ENDDO
      ENDIF
C
C  *** C36
C: CONSTANT BOTTOM AND SURFACE OBCS
C
      WRITE(2,999)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      DO M=1,5
        READ(1,90) TITLE(M)
        WRITE(2,90) TITLE(M)
      ENDDO
      NWQOBW=0
      M=0
      IF(NWQOBW_GL>0)THEN
        DO MM=1,NWQOBW_GL
          READ(1,*) I_TEMP,J_TEMP,(WQOBCW_GL(MM,1,NW),NW=1,NWQV)
          I_TEMP=XLOC(I_TEMP)
          J_TEMP=YLOC(J_TEMP)
          L=LIJ(I_TEMP,J_TEMP)
          IF (CELL_INSIDE_DOMAIN(L)) THEN
            NWQOBW=NWQOBW+1
            M=M+1
            IWQCBW(M) = I_TEMP
            JWQCBW(M) = J_TEMP
            WQOBCW(M,1,1:NWQV)=WQOBCW_GL(MM,1,1:NWQV)
          ! *** CONCENTRATION ASSIGNMENTS, BOTTOM
            IF(IWQCBW(M).EQ.ICBW(M).AND.JWQCBW(M).EQ.JCBW(M))THEN
              DO NW=1,NWQV
                NT=4+NTOX+NSED+NSND+NW
                CBW(M,1,NT)=WQOBCW(M,1,NW)
              ENDDO
            ENDIF
            WRITE(2,*) IWQCBW(M),JWQCBW(M),(WQOBCW(M,1,NW),NW=1,NWQV)
            WRITE(2,97) IWQCBW(M),JWQCBW(M),(WQOBCW(M,1,NW),NW=1,NWQV)
          ENDIF
        ENDDO
      ENDIF
C  *** C37
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      WRITE(2,999)
      DO M=1,5
        READ(1,90) TITLE(M)
        WRITE(2,90) TITLE(M)
      ENDDO
      NWQOBW=0
      M=0
      IF(NWQOBW_GL.GT.0)THEN
        DO MM=1,NWQOBW_GL
          READ(1,*) I_TEMP,J_TEMP,(WQOBCW_GL(MM,2,NW),NW=1,NWQV)
          I_TEMP=XLOC(I_TEMP)
          J_TEMP=YLOC(J_TEMP)
          L=LIJ(I_TEMP,J_TEMP)
          IF (CELL_INSIDE_DOMAIN(L)) THEN
            NWQOBW=NWQOBW+1
            M=M+1
            IWQCBW(M) = I_TEMP
            JWQCBW(M) = J_TEMP
            WQOBCW(M,2,1:NWQV)=WQOBCW_GL(MM,2,1:NWQV)
          ! *** CONCENTRATION ASSIGNMENTS, TOP
            IF(IWQCBW(M).EQ.ICBW(M).AND.JWQCBW(M).EQ.JCBW(M))THEN
              DO NW=1,NWQV
                NT=4+NTOX+NSED+NSND+NW
                CBW(M,2,NT)=WQOBCW(M,2,NW)
              ENDDO
            ENDIF
            WRITE(2,*) IWQCBW(M),JWQCBW(M),(WQOBCW(M,2,NW),NW=1,NWQV)
            WRITE(2,97) IWQCBW(M),JWQCBW(M),(WQOBCW(M,2,NW),NW=1,NWQV)
          ENDIF
        ENDDO
      ENDIF
C
C  *** C38
C EAST BOUNDARY
C
      WRITE(2,999)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      READ(1,90) TITLE(1)
      WRITE(2,90) TITLE(1)
      DO M=1,5
        READ(1,90) TITLE(M)
        WRITE(2,90) TITLE(M)
      ENDDO
      NWQOBE=0
      M=0
      IF(NWQOBE_GL.GT.0)THEN
        DO MM=1,NWQOBE_GL
          READ(1,*) I_TEMP,J_TEMP,(IWQOBE_GL(MM,NW),NW=1,NWQV)
          I_TEMP=XLOC(I_TEMP)
          J_TEMP=YLOC(J_TEMP)
          L=LIJ(I_TEMP,J_TEMP)
          IF (CELL_INSIDE_DOMAIN(L)) THEN
            NWQOBE=NWQOBE+1
            M=M+1
            IWQCBE(M) = I_TEMP
            JWQCBE(M) = J_TEMP
            IWQOBE(M,1:NWQV)=IWQOBE_GL(MM,1:NWQV)
          ! *** CONCENTRATION ASSIGNMENTS
            IF(IWQCBE(M).EQ.ICBE(M).AND.JWQCBE(M).EQ.JCBE(M))THEN
              DO NW=1,NWQV
                NT=4+NTOX+NSED+NSND+NW
                NCSERE(M,NT)=IWQOBE(M,NW)
              ENDDO
            ELSE
              STOP '  WQ: EAST OBC: MISMATCH BETWEEN NCBE & NWQOBE'
            ENDIF
            WRITE(2,*) IWQCBE(M),JWQCBE(M),(IWQOBE(M,NW),NW=1,NWQV)
            WRITE(2,969) IWQCBE(M),JWQCBE(M),(IWQOBE(M,NW),NW=1,NWQV)
          ENDIF
        ENDDO
      ENDIF
C
C  *** C39
C: CONSTANT BOTTOM AND SURFACE OBCS
C
      WRITE(2,999)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      DO M=1,5
        READ(1,90) TITLE(M)
        WRITE(2,90) TITLE(M)
      ENDDO
      NWQOBE=0
      M=0
      IF(NWQOBE_GL.GT.0)THEN
        DO MM=1,NWQOBE_GL
          READ(1,*) I_TEMP,J_TEMP,(WQOBCE_GL(MM,1,NW),NW=1,NWQV)
          I_TEMP=XLOC(I_TEMP)
          J_TEMP=YLOC(J_TEMP)
          L=LIJ(I_TEMP,J_TEMP)
          IF (CELL_INSIDE_DOMAIN(L)) THEN
            NWQOBE=NWQOBE+1
            M=M+1
            IWQCBE(M) = I_TEMP
            JWQCBE(M) = J_TEMP
            WQOBCE(M,1,1:NWQV)=WQOBCE_GL(MM,1,1:NWQV)
          ! *** CONCENTRATION ASSIGNMENTS, BOTTOM
            IF(IWQCBE(M).EQ.ICBE(M).AND.JWQCBE(M).EQ.JCBE(M))THEN
              DO NW=1,NWQV
                NT=4+NTOX+NSED+NSND+NW
                CBE(M,1,NT)=WQOBCE(M,1,NW)
              ENDDO
            ENDIF
            WRITE(2,*) IWQCBE(M),JWQCBE(M),(WQOBCE(M,1,NW),NW=1,NWQV)
            WRITE(2,97) IWQCBE(M),JWQCBE(M),(WQOBCE(M,1,NW),NW=1,NWQV)
          ENDIF
        ENDDO
      ENDIF
      WRITE(2,999)
C  *** C40
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      DO M=1,5
        READ(1,90) TITLE(M)
        WRITE(2,90) TITLE(M)
      ENDDO
      NWQOBE=0
      M=0
      IF(NWQOBE_GL.GT.0)THEN
        DO MM=1,NWQOBE_GL
          READ(1,*) I_TEMP,J_TEMP,(WQOBCE_GL(MM,2,NW),NW=1,NWQV)
          I_TEMP=XLOC(I_TEMP)
          J_TEMP=YLOC(J_TEMP)
          L=LIJ(I_TEMP,J_TEMP)
          IF (CELL_INSIDE_DOMAIN(L)) THEN
            NWQOBE=NWQOBE+1
            M=M+1
            IWQCBE(M) = I_TEMP
            JWQCBE(M) = J_TEMP
            WQOBCE(M,2,1:NWQV)=WQOBCE_GL(MM,2,1:NWQV)
          ! *** CONCENTRATION ASSIGNMENTS, TOP
            IF(IWQCBE(M).EQ.ICBE(M).AND.JWQCBE(M).EQ.JCBE(M))THEN
              DO NW=1,NWQV
                NT=4+NTOX+NSED+NSND+NW
                CBE(M,2,NT)=WQOBCE(M,2,NW)
              ENDDO
            ENDIF
            WRITE(2,*) IWQCBE(M),JWQCBE(M),(WQOBCE(M,2,NW),NW=1,NWQV)
            WRITE(2,97) IWQCBE(M),JWQCBE(M),(WQOBCE(M,2,NW),NW=1,NWQV)
          ENDIF
        ENDDO
      ENDIF
C
C  *** C41
C NORTH BDRY
C       READ(1,90) TITLE(M)
C
      WRITE(2,999)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      READ(1,90) TITLE(1)
      WRITE(2,90) TITLE(1)
      DO M=1,5
        READ(1,90) TITLE(M)
        WRITE(2,90) TITLE(M)
      ENDDO
      NWQOBN=0
      M=0
      IF(NWQOBN_GL.GT.0)THEN
        DO MM=1,NWQOBN_GL
          READ(1,*) I_TEMP,J_TEMP,(IWQOBN_GL(MM,NW),NW=1,NWQV)
          I_TEMP=XLOC(I_TEMP)
          J_TEMP=YLOC(J_TEMP)
          L=LIJ(I_TEMP,J_TEMP)
          IF (CELL_INSIDE_DOMAIN(L)) THEN
            NWQOBN=NWQOBN+1
            M=M+1
            IWQCBN(M) = I_TEMP
            JWQCBN(M) = J_TEMP
            IWQOBN(M,1:NWQV)=IWQOBN_GL(MM,1:NWQV)
          ! *** CONCENTRATION ASSIGNMENTS
            IF(IWQCBN(M).EQ.ICBN(M).AND.JWQCBN(M).EQ.JCBN(M))THEN
              DO NW=1,NWQV
                NT=4+NTOX+NSED+NSND+NW
                NCSERN(M,NT)=IWQOBN(M,NW)
              ENDDO
            ELSE
              STOP '  WQ: NORTH OBC: MISMATCH BETWEEN NCBN & NWQOBN'
            ENDIF
            WRITE(2,*) IWQCBN(M),JWQCBN(M),(IWQOBN(M,NW),NW=1,NWQV)
            WRITE(2,969) IWQCBN(M),JWQCBN(M),(IWQOBN(M,NW),NW=1,NWQV)
          ENDIF
        ENDDO
      ENDIF
C
C  *** C42
C: CONSTANT BOTTOM AND SURFACE OBCS
C
      WRITE(2,999)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      DO M=1,5
        READ(1,90) TITLE(M)
        WRITE(2,90) TITLE(M)
      ENDDO
      NWQOBN=0
      M=0
      IF(NWQOBN_GL.GT.0)THEN
        DO MM=1,NWQOBN_GL
          READ(1,*) I_TEMP,J_TEMP,(WQOBCN_GL(MM,1,NW),NW=1,NWQV)
          I_TEMP=XLOC(I_TEMP)
          J_TEMP=YLOC(J_TEMP)
          L=LIJ(I_TEMP,J_TEMP)
          IF (CELL_INSIDE_DOMAIN(L)) THEN
            NWQOBN=NWQOBN+1
            M=M+1
            IWQCBN(M) = I_TEMP
            JWQCBN(M) = J_TEMP
            WQOBCN(M,1,1:NWQV)=WQOBCN_GL(MM,1,1:NWQV)
          ! *** CONCENTRATION ASSIGNMENTS, BOTTOM
            IF(IWQCBN(M).EQ.ICBN(M).AND.JWQCBN(M).EQ.JCBN(M))THEN
              DO NW=1,NWQV
                NT=4+NTOX+NSED+NSND+NW
                CBN(M,1,NT)=WQOBCN(M,1,NW)
              ENDDO
            ENDIF
            WRITE(2,*) IWQCBN(M),JWQCBN(M),(WQOBCN(M,1,NW),NW=1,NWQV)
            WRITE(2,97) IWQCBN(M),JWQCBN(M),(WQOBCN(M,1,NW),NW=1,NWQV)
          ENDIF
        ENDDO
      ENDIF
      WRITE(2,999)
C  *** C43
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      DO M=1,5
        READ(1,90) TITLE(M)
        WRITE(2,90) TITLE(M)
      ENDDO
      NWQOBN=0
      M=0
      IF(NWQOBN_GL.GT.0)THEN
        DO MM=1,NWQOBN_GL
          READ(1,*) I_TEMP,J_TEMP,(WQOBCN_GL(MM,2,NW),NW=1,NWQV)
          I_TEMP=XLOC(I_TEMP)
          J_TEMP=YLOC(J_TEMP)
          L=LIJ(I_TEMP,J_TEMP)
          IF (CELL_INSIDE_DOMAIN(L)) THEN
            NWQOBN=NWQOBN+1
            M=M+1
            IWQCBN(M) = I_TEMP
            JWQCBN(M) = J_TEMP
            WQOBCN(M,2,1:NWQV)=WQOBCN_GL(MM,2,1:NWQV)
          ! *** CONCENTRATION ASSIGNMENTS, TOP
            IF(IWQCBN(M).EQ.ICBN(M).AND.JWQCBN(M).EQ.JCBN(M))THEN
              DO NW=1,NWQV
                NT=4+NTOX+NSED+NSND+NW
                CBN(M,2,NT)=WQOBCN(M,2,NW)
              ENDDO
            ENDIF
            WRITE(2,*) IWQCBN(M),JWQCBN(M),(WQOBCN(M,2,NW),NW=1,NWQV)
            WRITE(2,97) IWQCBN(M),JWQCBN(M),(WQOBCN(M,2,NW),NW=1,NWQV)
          ENDIF
        ENDDO
      ENDIF
C
C  *** C44
C SPATIALLY/TEMPORALLY CONSTANT INITIAL CONDITIONS: WQCHLX=1/WQCHLX
C READ DATA POINTS & DO INTERNAL INTERPOLATION?
C
      WRITE(2,999)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      READ(1,90) TITLE(1)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0)THEN
        DO M=1,3
          READ(1,999)
        ENDDO
      ENDIF
      READ(1,*) (WQV(1,1,NW), NW=1,6)
      WRITE(2,*) (WQV(1,1,NW), NW=1,6)
      READ(1,*) (WQV(1,1,NW), NW=7,13)
      WRITE(2,*) (WQV(1,1,NW), NW=7,13)
      IF(IDNOTRVA/=0)THEN
        READ(1,*) (WQV(1,1,NW), NW=14,NWQV),WQV(1,1,IDNOTRVA),WQMCMIN
        WRITE(2,*) (WQV(1,1,NW), NW=14,NWQV),WQV(1,1,IDNOTRVA),WQMCMIN
      ELSE
        READ(1,*) (WQV(1,1,NW), NW=14,NWQV)
        WRITE(2,*) (WQV(1,1,NW), NW=14,NWQV)
      ENDIF
      IF(IWQICI.EQ.0)THEN  !SCJ only needed if ICs are constant and specified on C44
        WRITE(2,999)
        WRITE(2,90) TITLE(1)
        WRITE(2,21)' : (BC, BD, BG)         = ', (WQV(1,1,NW),NW=1,3)
        WRITE(2,21)' : (RPOC, LPOC, DOC)    = ', (WQV(1,1,NW),NW=4,6)
        WRITE(2,21)' : (RPOP,LPOP,DOP,PO4T) = ', (WQV(1,1,NW),NW=7,10)
        WRITE(2,21)' : (RPON, LPON, DON)    = ', (WQV(1,1,NW),NW=11,13)
        WRITE(2,21)' : (NH4, NO3)           = ', (WQV(1,1,NW),NW=14,15)
        WRITE(2,21)' : (SU, SA, COD, DO)    = ', (WQV(1,1,NW),NW=16,19)
        WRITE(2,981)' : (TAM, FCB, CO2, MAC)= ',
     &		                (WQV(1,1,NW),NW=20,NWQV),WQV(1,1,IDNOTRVA)			!VB CHANGED NWQV+1 TO NWQV 	!VB added CO2
     &																			!ADDED WQV(1,1,IDNOTRVA)
        WQCHL(1,1) = WQV(1,1,1)*WQCHLC + WQV(1,1,2)*WQCHLD
     &      + WQV(1,1,3)*WQCHLG
        IF(IWQSRP.EQ.1)THEN
          O2WQ_ = MAX(WQV(1,1,19), 0.0)
          WQTAMD = MIN( WQTAMDMX*EXP(-WQKDOTAM*O2WQ_), WQV(1,1,20) )
          WQTAMP(1,1) = WQV(1,1,20) - WQTAMD
          WQPO4D(1,1) = WQV(1,1,10) / (1.0 + WQKPO4P*WQTAMP(1,1))
          WQSAD(1,1)  = WQV(1,1,17) / (1.0 + WQKSAP*WQTAMP(1,1))
        ELSEIF(IWQSRP.EQ.2)THEN
          WQPO4D(1,1) = WQV(1,1,10) / (1.0 + WQKPO4P*SEDT(1,1))
          WQSAD(1,1)  = WQV(1,1,17) / (1.0 + WQKSAP*SEDT(1,1))
        ELSE
          WQPO4D(1,1) = WQV(1,1,10)
          WQSAD(1,1)  = WQV(1,1,17)
        ENDIF
        DO NW=1,NWQV
          TVARWQ=WQV(1,1,NW)
          DO K=1,KC
            WQV(LC,K,NW) = TVARWQ
            WQV(1 ,K,NW) = TVARWQ
            WQVO(LC,K,NW) = TVARWQ
            WQVO(1 ,K,NW) = TVARWQ
          ENDDO
        ENDDO
        DO NW=1,NWQV
          TVARWQ=WQV(1,1,NW)
          DO ND=1,NDMWQ
            LF=2+(ND-1)*LDMWQ
            LL=LF+LDM-1
            DO K=1,KC
              DO L=LF,LL
                WQV(L,K,NW) = TVARWQ
              ENDDO
            ENDDO
          ENDDO
        ENDDO
        DO NW=1,NWQV !These are required to initialize WQ variables along boundaries so that concentration discontinuities can be smoothed.
          M=4+NTOX+NSED+NSND+NW
          DO K=1,KC  !If there are data on EFDC.INP C47, C52, C57, C62 and WQ is active, these need to be initialized
            DO LL=1,NWQOBS
              CLOS(LL,K,M)=WQV(1,1,NW)
              NLOS(LL,K,M)=N
            ENDDO
            DO LL=1,NWQOBW
              CLOW(LL,K,M)=WQV(1,1,NW)
              NLOW(LL,K,M)=N
            ENDDO
            DO LL=1,NWQOBE
              CLOE(LL,K,M)=WQV(1,1,NW)
              NLOE(LL,K,M)=N
            ENDDO
            DO LL=1,NWQOBN
              CLON(LL,K,M)=WQV(1,1,NW)
              NLON(LL,K,M)=N
            ENDDO
          ENDDO
        ENDDO
        XWQCHL = WQCHL(1,1)
        XWQTAMP = WQTAMP(1,1)
        XWQPO4D = WQPO4D(1,1)
        XWQSAD = WQSAD(1,1)
        DO K=1,KC
          WQCHL(LC,K) = XWQCHL
          WQTAMP(LC,K) = XWQTAMP
          WQPO4D(LC,K) = XWQPO4D
          WQSAD(LC,K) = XWQSAD
          WQCHL(1,K) = XWQCHL
          WQTAMP(1,K) = XWQTAMP
          WQPO4D(1,K) = XWQPO4D
          WQSAD(1,K) = XWQSAD
        ENDDO
        DO ND=1,NDMWQ
          LF=2+(ND-1)*LDMWQ
          LL=LF+LDM-1
          DO K=1,KC
            DO L=LF,LL
              WQCHL(L,K) = XWQCHL
              WQTAMP(L,K) = XWQTAMP
              WQPO4D(L,K) = XWQPO4D
              WQSAD(L,K) = XWQSAD
            ENDDO
          ENDDO
        ENDDO
      ENDIF
C
C READ CELL MAPPING FILE 'MACALGMP.INP' AND SET INITIAL CONDITION.
C MACROALGAE WILL RESIDE ONLY IN THE BOTTOM LAYER.
C 09/02/99 M.MORTON: ADDED KMV VELOCITY HALF-SATURATION, KBP DENSITY
C   HALF SATURATION, AND PSHADE FACTOR TO BETTER CONTROL MACROALGAE
C   GROWTH.
C SCJ - Macroalgae allowed in all layers (not just bottom). MACALGMP.INP only read if there is growth-limitation due to velocities or biomass.
C
      IF(IWQVLIM>0)THEN
        DO L=1,LCM
          SMAC(L)=0.0
          PSHADE(L)=1.0
          WQKMV(L)=0.25
          WQKMVMIN(L)=0.15
          WQKBP(L)=6.5
          WQKMVA(L)=1.0
          WQKMVB(L)=12.0
          WQKMVC(L)=0.3
          WQKMVD(L)=0.25
          WQKMVE(L)=2.0
        ENDDO
        WRITE(2,9003)
 9003 FORMAT(/,' MACALGMP.INP - MACROALGAE MAP FILE',/,
     &    ' PSHADE = SHADE FACTOR FOR TREE CANOPY (1.0=NO CANOPY)',/,
     &    ' KMV    = MACROALGAE HALF-SATURATION VELOCITY (M/SEC)',/,
     &    ' KMVMIN = MACROALGAE VELOCITY LIMITATION MINIMUM (M/SEC)',/,
     &    ' KBP    = MACROALGAE HALF-SATURATION DENSITY (G C/M2)',/,
     &    ' KMVA   = MACROALGAE VEL. LIMIT LOGISTIC FUNC. PARAM. A',/,
     &    ' KMVB   = MACROALGAE VEL. LIMIT LOGISTIC FUNC. PARAM. B',/,
     &    ' KMVC   = MACROALGAE VEL. LIMIT LOGISTIC FUNC. PARAM. C',/,
     &    ' KMVD   = MACROALGAE VEL. LIMIT LOGISTIC FUNC. PARAM. D',/,
     &    ' KMVE   = MACROALGAE VEL. LIMIT LOGISTIC FUNC. PARAM. E',/,
     &    '   I   J   L PSHADE    KMV KMVMIN    KBP   KMVA   KMVB',
     &    '   KMVC   KMVD   KMVE')
        PRINT *,'WQ: MACALGMP.INP'
        OPEN(3,FILE='MACALGMP.INP',STATUS='UNKNOWN')
        CALL SKIPCOMM(3, CCMRM)
 9001   READ(3,*,END=9002) I, J, XMRM1, XMRM2, XMRM3, XMRM4,
     &      XMRMA, XMRMB, XMRMC, XMRMD, XMRME
        II = XLOC(I)
        JJ = YLOC(J)
        LL=LIJ(II,JJ)
        IF(II .LE. 0) GOTO 9002
        IF (CELL_INSIDE_DOMAIN(LL)) THEN
          IF(IJCT(II,JJ).LT.1 .OR. IJCT(II,JJ).GT.8)THEN
            PRINT*, 'I, J, IJCT(I,J) = ', II,JJ,IJCT(II,JJ)
            STOP 'ERROR!! INVALID (I,J) IN FILE MACALGMP.INP'
          ENDIF
          SMAC(LL)=1.0
          PSHADE(LL)=XMRM1    ! *** PMC-This overwrites the Heat Module Shading
          WQKMV(LL)=XMRM2
          WQKMVMIN(LL)=XMRM3
          WQKBP(LL)=XMRM4
          WQKMVA(LL)=XMRMA
          WQKMVB(LL)=XMRMB
          WQKMVC(LL)=XMRMC
          WQKMVD(LL)=XMRMD
          WQKMVE(LL)=XMRME
          WRITE(2,9004) II, JJ, LL, PSHADE(LL), WQKMV(LL), WQKMVMIN(LL),
     &      WQKBP(LL), WQKMVA(LL), WQKMVB(LL), WQKMVC(LL), WQKMVD(LL),
     &      WQKMVE(LL)
 9004 FORMAT(' ',I3,' ',I3,' ',I5, 9F7.3)
          GOTO 9001
 9002     CLOSE(3)
        ENDIF
      ENDIF
      !INITIALIZE MACROALGAE CONCENTRATIONS TO VEGE.INP PRESCRIBED CONDITIONS
      IF(LOGMAC)THEN
        DO M=1,MCOUNT
          L=IJLMAC(M,3)
          DO K=1,KC !Make sure it is in more than the bottom layer
            !FRACLAY=FLOAT(K)/FLOAT(KC)  !top of the layer (normalized)
            FRACLAY=SUM(DZC(1:K)) !top of the layer as a fraction of HP
            !FRACLAYKM1=FLOAT(K-1)/FLOAT(KC)  !top of the layer (normalized)
            FRACLAYKM1=SUM(DZC(1:K-1)) !top of the layer as a fraction of HP
            FHLAYC=FRACLAY*HP(L)  !z location of the normalized top of the sigma layer of the cell
            FHLAYCKM1=FRACLAYKM1*HP(L)  !z location of the normalized top of the sigma layer of the cell
            WQVMTMP=(0.25*PI*(MACDIAM(L)*MACDIAM(L)))*RDLPSQ(MVEGL(L))*1027.0*0.14*0.22*1000.0 ! Total grams-Carbon of WQV macroalgae per volume water
! WQV(C/m3) = [cross-sectional area of single macroalgae branch] x [number of branches per unit area] x [density of seawater in kg x percentage of algae dry weight to total undried weight x percentage of carbon to total dry weight] x [1000 g per kg to convert from kg to g]
            !PRINT *, 'WQVMTMP, RDLPSQ(M)) = ', WQVMTMP, RDLPSQ(MVEGL(L))
            !...multiplied by neutral density and 90%/10% water:biomass content and 36% carbon content by dry weight
            IF(ZMAXMAC(L)>=FHLAYCKM1)THEN
              IF(ZMAXMAC(L)<FHLAYC)WQV(L,K,IDNOTRVA)=WQVMTMP*(ZMAXMAC(L)-FHLAYCKM1)/(FHLAYC-FHLAYCKM1) ! scaled by partial fill ratio
              IF(ZMAXMAC(L)>=FHLAYC)WQV(L,K,IDNOTRVA)=WQVMTMP
            ELSE
              WQV(L,K,IDNOTRVA)=0.0
            ENDIF
            IF(ZMINMAC(L)>FHLAYC)THEN    ! Algae not in this layer
              WQV(L,K,IDNOTRVA)=0.0
            ELSE
              IF(ZMINMAC(L)>=FHLAYCKM1)THEN    !This single layer completely filled with macroalgae
                IF(ZMAXMAC(L)>=FHLAYC)THEN
                  WQV(L,K,IDNOTRVA)=WQVMTMP*(FHLAYC-ZMINMAC(L))/(FHLAYC-FHLAYCKM1)
                ELSEIF(ZMAXMAC(L)<FHLAYC)THEN
                  WQV(L,K,IDNOTRVA)=WQVMTMP*(ZMAXMAC(L)-ZMINMAC(L))/(FHLAYC-FHLAYCKM1)
                ENDIF
              ENDIF
            ENDIF
            WQVO(L,K,IDNOTRVA)=WQV(L,K,IDNOTRVA)
          ENDDO
        ENDDO
      ENDIF
C
C
C  *** C45
C SPATIALLY/TEMPORALLY CONSTANT ALGAL GROWTH, RESPIRATION & PREDATION RA
C
      WRITE(2,999)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      READ(1,90) TITLE(1)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      READ(1,*) WQPMC(1),WQPMD(1),WQPMG(1),WQPMM(1),WQBMRC(1),
     &    WQBMRD(1),WQBMRG(1),WQBMRM(1),WQPRRC(1),WQPRRD(1),
     &    WQPRRG(1),WQPRRM(1),WQKEB(1)
      WRITE(2,*) WQPMC(1),WQPMD(1),WQPMG(1),WQPMM(1),WQBMRC(1),
     &    WQBMRD(1),WQBMRG(1),WQBMRM(1),WQPRRC(1),WQPRRD(1),
     &    WQPRRG(1),WQPRRM(1),WQKEB(1)
      IF(IWQAGR.NE.1)THEN
        WRITE(2,999)
        WRITE(2,90) TITLE(1)
      WRITE(2,80)'* ALGAL GROWTH RATE (/DAY)                        '
      WRITE(2,21)' : (PMC, PMD, PMG)       = ', WQPMC(1),WQPMD(1),
     &    WQPMG(1)
      WRITE(2,80)'* ALGAL BASAL METABOLISM RATE (/DAY)              '
      WRITE(2,21)' : (BMRC, BMRD, BMRG)    = ', WQBMRC(1),WQBMRD(1),
     &    WQBMRG(1)
      WRITE(2,80)'* ALGAL PREDATION RATE (/DAY)                     '
      WRITE(2,21)' : (PRRC, PRRD, PRRG)    = ', WQPRRC(1),WQPRRD(1),
     &    WQPRRG(1)
        WRITE(2,82)
     &      '* BASE LIGHT EXTINCTION COEFFICIENT (/M)   = ',WQKEB(1)
        DO I=2,IWQZ
          WQPMC(I)=WQPMC(1)
          WQPMD(I)=WQPMD(1)
          WQPMG(I)=WQPMG(1)
          WQPMM(I)=WQPMM(1)
          WQBMRC(I)=WQBMRC(1)
          WQBMRD(I)=WQBMRD(1)
          WQBMRG(I)=WQBMRG(1)
          WQBMRM(I)=WQBMRM(1)
          WQPRRC(I)=WQPRRC(1)
          WQPRRD(I)=WQPRRD(1)
          WQPRRG(I)=WQPRRG(1)
          WQPRRM(I)=WQPRRM(1)
          WQKEB(I)=WQKEB(1)
        ENDDO
      ENDIF
C
C *** C46 SPATIALLY/TEMPORALLY CONSTANT SETTLING VELOCITIES AND REAERATION FACTO
C
      WRITE(2,999)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      READ(1,90) TITLE(1)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      READ(1,*)WQWSC(1),WQWSD(1),WQWSG(1),WQWSRP(1),WQWSLP(1),WQWSS(1),
     &    WQWSM, REAC(1)
      WRITE(2,*)WQWSC(1),WQWSD(1),WQWSG(1),WQWSRP(1),WQWSLP(1),
     &    WQWSS(1),WQWSM, REAC(1)
      IF(IWQSTL.NE.1)THEN
        WRITE(2,999)
        WRITE(2,90) TITLE(1)
      WRITE(2,80)'* ALGAL SETTLING RATE (M/DAY)                     '
      WRITE(2,21)' : (WSC, WSD, WSG)       = ', WQWSC(1),WQWSD(1),WQWSG(1)
      WRITE(2,80)'* POM SETTLING RATE (M/DAY)                       '
      WRITE(2,21)' : (WSRP, WSLP)          = ', WQWSRP(1),WQWSLP(1)
      WRITE(2,80)'* SETTLING RATE OF PARTICULATE METAL (M/DAY)      '
      WRITE(2,21)' : (WSS)                 = ', WQWSS(1)
        DO I=2,IWQZ
          WQWSC(I)=WQWSC(1)
          WQWSD(I)=WQWSD(1)
          WQWSG(I)=WQWSG(1)
          WQWSRP(I)=WQWSRP(1)
          WQWSLP(I)=WQWSLP(1)
          WQWSS(I)=WQWSS(1)
          REAC(I)=REAC(1)
        ENDDO
      ENDIF
C
C *** C47 SPATIALLY/TEMPORALLY CONSTANT BENTHIC FLUXES
C
      WRITE(2,999)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      READ(1,90) TITLE(1)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      READ(1,*) WQBFPO4D(1),WQBFNH4(1),WQBFNO3(1),WQBFSAD(1),WQBFCOD(1),WQBFO2(1)
      WRITE(2,*) WQBFPO4D(1),WQBFNH4(1),WQBFNO3(1),WQBFSAD(1),WQBFCOD(1),WQBFO2(1)
      IF(IWQBEN.EQ.0)THEN
        WRITE(2,999)
        WRITE(2,90) TITLE(1)
      WRITE(2,21)' : (PO4D, NH4, NO3)     = ',WQBFPO4D(1),WQBFNH4(1),WQBFNO3(1)
      WRITE(2,21)' : (SAD, COD, DO)       = ',WQBFSAD(1),WQBFCOD(1),WQBFO2(1)
        DO L=2,LA
          WQBFPO4D(L)=WQBFPO4D(1)
          WQBFNH4(L)=WQBFNH4(1)
          WQBFNO3(L)=WQBFNO3(1)
          WQBFSAD(L)=WQBFSAD(1)
          WQBFCOD(L)=WQBFCOD(1)
          WQBFO2(L)=WQBFO2(1)
        ENDDO
      ENDIF
C
C *** C48
C
C *** TEMPORALLY-CONSTANT VALUES FOR POINT SOURCE CONCENTRATIONS IN MG/L
C ***  EXCEPT XPSQ (M^3/S), TAM (KMOL/D), FCB (MPN/L).
C
      WRITE(2,999)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      READ(1,90) TITLE(1)
      WRITE(2,999)
      WRITE(2,90) TITLE(1)
C *** C48
C     This is where constant point source loads are specified.
C     The code was updated to be on a layer-wise basis for macroalgae "irrigation."
C     Original code ignored K, but it is used for irrigating macroalgae.
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      IF(ISSKIP .EQ. 0) READ(1,999)
      READ(1,*) IWQPS_GL,NPSTMSR
      WRITE(2,*) IWQPS_GL,NPSTMSR
      WRITE(2,23)'* NUMBER OF CELLS FOR POINT SOURCE INPUT  = ',IWQPS_GL
      WRITE(2,23)'* NUMBER WITH VARIABLE POINT SOURCE INPUT = ',NPSTMSR
      IF(IWQPS_GL.GT.NWQPS_GL) STOP 'ERROR!! IWQPS SHOULD BE <= NWQPS'
      DO M=1,3
        READ(1,90) TITLE(M)
        WRITE(2,90) TITLE(M)
      ENDDO
      M=0
      IWQPS=0
      WQPSLs: DO MM=1,IWQPS_GL
        READ(1,*) I,J,K,ITMP,XPSQ,(XPSL(NW),NW=1,6)  !GR changed ITMP to be indexed value CSER in CWQSRxx.INP, K is sigma layer for inflow corresponding to QSER.INP. XPSQ is not used.
        WRITE(2,*) I,J,K,ITMP,XPSQ,(XPSL(NW),NW=1,6) !SCJ ITMP (formerly KSR) is the flow series index in QSER.INP, but C32, C35, C38, and C41 already do this for BCs (not irrigation)
        READ(1,*) (XPSL(NW),NW=7,13)
        WRITE(2,*) (XPSL(NW),NW=7,13)
        READ(1,*) (XPSL(NW),NW=14,NWQV)
        WRITE(2,*) (XPSL(NW),NW=14,NWQV)
        WRITE(2,294) I,J,K,ITMP,XPSQ,(XPSL(NW),NW=1,NWQV)
        II = XLOC(I)
        JJ = YLOC(J)
        L = LIJ(II,JJ)
        !*********************************NOTE
        !If we want to use C48 to specify irrigation, we ned to make sure CSER is only applied at layer K in C48.
        !Presently it is assigned to all layers in CALCSER.
        IF(CELL_INSIDE_DOMAIN(L))THEN
          M=M+1
          IWQPS = IWQPS + 1
          IF(MOD(M,NQSIJ+1)==0)M=1 !Reset M to allow distinct layers to be different (not assigned to the entire water column)
          IF(IJCT(II,JJ)<1 .OR. IJCT(II,JJ)>8)THEN
            WRITE(*,911) I,J
            STOP 'ERROR!! INVALID (I,J) IN FILE WQ3DWC.INP FOR PSL'
          ENDIF
        ! *** HANDLE CONCENTRATION-BASED POINT SOURCE
          IF(IWQPSL.EQ.2.OR.IWQPSL.EQ.3)THEN !GR Feb 2019 - added L,K indexed time series volume capabilities to WQV using C48 inupts and IWQPSL=3
            IF(L.NE.LQS(M))THEN
              STOP ' MISMATCH NQSIJ BETWEEN EFDC.INP & WQ3DWC.INP'
            ENDIF
          ! *** ASSIGN GLOBAL CONCENTRATION TIME SERIES INDEX
            DO NW=1,NWQV
              N1=4+NTOX+NSED+NSND+NW
! NCSERQ is referenced in CALFQC during standard volume source/sinks *** SCJ set WQ time series
              NCSERQ(M,N1)=ITMP  ! *** For macroalgae irrigation, use the flow time series in QSER.INP. All WQ variables use the same flow time series.
            ! *** CONVERT FROM Kmol TO moles
            !WQWPSLC(M,20) = XPSL(20) * CONV1   PMC?
            ! *** CONVERT FROM MPN/L TO MPN/DAY
            !WQWPSLC(M,NWQV) = XPSL(NWQV) * WQTT * CONV1   PMC?
              IF(IWQPSL==2)THEN
                CQS(1:KC,M,N1)=XPSL(NW)  !Existing EFDC code where K was ignored in C48
              ELSE !K is read from C48
                CQS(K   ,M,N1)=XPSL(NW)!Set constant point load contribution. Set to XPSL values in WQ3DWC.INP to 0 if you want to only use time series in CWQSRxx.INP
              ENDIF
            ENDDO
          ENDIF

C
C JMH MODIFIED 5/18/00 TO ALLOW KCPSL(M) TO BE SET TO ZERO FOR UNIFORM D
C OF LOAD IN HORIZONTAL CELL STACK OVER ALL LAYERS
C
          DO KK=1,KC
            IWQPSC(L,KK)=M
            IWQPSV(L,KK)=ITMP
          ENDDO
          ICPSL(M)=II
          JCPSL(M)=JJ
          KCPSL(M)=K
          MVPSL(M)=ITMP
C        WQPSQC(M)=XPSQ1  ! *** PMC-NOT USED

        ! *** CONVERT FROM CONC (mg/l) AND Q (m3/s) TO MASS (G/DAY)
        ! *** CONV1=1.0  AND   CONV2=8.64E4.
          XPSQ1=-XPSQ !need this variable for two things depending on whether IWQPSL=1 or 2
          WQTT = -XPSQ1*CONV2
          DO NW=1,19
            WQWPSLC(M,NW) = XPSL(NW) * CONV1 * WQTT
          ENDDO
        ! *** CONVERT FROM Kmol TO moles
          WQWPSLC(M,20) = XPSL(20) * CONV1
        ! *** CONVERT FROM MPN/L TO MPN/DAY
          WQWPSLC(M,21) = XPSL(21) * WQTT * CONV1
        ! *** CO2
          WQWPSLC(M,22) = XPSL(22) * CONV1 * WQTT
	    IF(IWQPSL==2)THEN
            FORALL(KK=1:KC)WQWPSL(L,KK,1:NWQV) !Be careful, this assumes that the layer was not specified in C48, but we may be changing this.
     &         =WQWPSL(L,KK,1:NWQV)+WQWPSLC(M,1:NWQV)*DZC(KK)/HP(L) !SCJ changed to have constant point source loads added correctly
          ELSE
             WQWPSL(L,KCPSL(M),1:NWQV)=WQWPSL(L,KCPSL(M),1:NWQV)+WQWPSLC(M,1:NWQV) !This is appropriate when K in C48 is used as the layer number
          ENDIF
        ENDIF
      ENDDO WQPSLs
      IF(NQSIJ_GL/=0.AND.IWQPS_GL/=0)THEN
        IF(IWQPS/=NQSIJ)THEN !Allows distinct layers to be different, not just assigned by column.
          PRINT*,'IWQPS not equal to NQSIJ'
          PRINT*,ANS(PARTID2),IWQPS,NQSIJ
          PAUSE
          STOP
        ENDIF
      ENDIF
  911 FORMAT(/,' I,J = ', 2I5)
C
C *** C49
      !  1) CHC - cyanobacteria
      !  2) CHD - diatom algae
      !  3) CHG - green algae
      !  4) ROC - refractory particulate organic carbon
      !  5) LOC - labile particulate organic carbon
      !  6) DOC - dissolved organic carbon
      !  7) ROP - refractory particulate organic phosphorus
      !  8) LOP - labile particulate organic phosphorus
      !  9) DOP - dissolved organic phosphorus
      ! 10) P4D - total phosphate
      ! 11) RON - refractory particulate organic nitrogen 22) macroalgae
      ! 12) LON - labile particulate organic nitrogen
      ! 13) DON - dissolved organic nitrogen
      ! 14) NHX - ammonia nitrogen
      ! 15) NOX - nitrate nitrogen
      ! 16) SUU - particulate biogenic silica
      ! 17) SAA - dissolved available silica
      ! 18) COD - chemical oxygen demand
      ! 19) DOX - dissolved oxygen
      ! 20) TAM - total active metal
      ! 21) FCB - fecal coliform bacteria
      ! 22) CO2 - dissolved CO2
C
C *** SPATIALLY/TEMPORALLY-CONSTANT VALUES FOR NON-POINT SOURCE INPUT
C ***   CONSTITUENT UNITS OF G/M2/DAY EXCEPT FCB WHICH IS MPN/M2/DAY.
C
      WRITE(2,999)
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      write(*,*) 'BEGIN DEBUG'
      READ(1,90) TITLE(1)
      write(*,90) TITLE(1)
      DO M=1,3
        READ(1,999)
      ENDDO
      READ(1,*) XDSQ,(XDSL(NW),NW=1,6)
      WRITE(2,*) XDSQ,(XDSL(NW),NW=1,6)
      READ(1,*) (XDSL(NW),NW=7,13)
      WRITE(2,*) (XDSL(NW),NW=7,13)
      READ(1,*) (XDSL(NW),NW=14,NWQV)
      WRITE(2,*) (XDSL(NW),NW=14,NWQV)
      IF(IWQNPL.NE.1)THEN  !SPATIALLY/TEMPORALLY CONSTANT NPS INPUT
        WRITE(2,999)
        WRITE(2,90) TITLE(1)
        WRITE(2, 21)' : (DSQ, CHC, CHD, CHG)  = ',XDSQ,(XDSL(NW),NW=1,3)
        WRITE(2, 21)' : (ROC, LOC, DOC)       = ',(XDSL(NW),NW=4,6)
        WRITE(2, 21)' : (ROP, LOP, DOP, P4D)  = ',(XDSL(NW),NW=7,10)
        WRITE(2, 21)' : (RON, LON, DON)       = ',(XDSL(NW),NW=11,13)
        WRITE(2, 21)' : (NHX, NOX)            = ',(XDSL(NW),NW=14,15)
        WRITE(2, 21)' : (SUU, SAA, COD, DOX)  = ',(XDSL(NW),NW=16,19)
        WRITE(2,981)' : (TAM, FCB, CO2)  = ',(XDSL(NW),NW=20,NWQV)    !VB ADDED co2
C PMC        WQDSQ(1,1) = XDSQ
C PMC        DO NW=1,18
C PMC          WQWDSL(1,1,NW) = XDSL(NW) * CONV1  ! CONVERT FROM KG/DAY TO G/DAY
C PMC        ENDDO

        ! *** CONVERT FROM Kmol TO moles
        WQTT = XDSQ*CONV2  !converts from seconds to days

        ! *** COMPUTE MASS LOADINGS, G/DAY, MOLES/DAY & MPN/DAY
        DO L=2,LA
          DO NW=1,NWQV
C
C M. MORTON MODIFIED THE LINE BELOW SO THAT CONSTANT ATMOSPHERIC DEPOSIT
C CAN BE ADDED VIA THIS ROUTINE INSTEAD OF CONSTANT NPS INPUT WHICH THE
C ORIGINAL CODE CALLED FOR AND WHICH WAS NOT PARTICULARLY USEFUL.
C INPUT DATA (XDSL) ARE IN G/M2/DAY AND ARE MULTIPLIED BY THE CELL SURFA
C AREA (DXYP) TO GET G/DAY.  ATMOSPHERIC DEPOSITION ONLY ENTERS THRU SUR
C LAYER (KC):
C
            WQWDSL(L,KC,NW) = XDSL(NW) * DXYP(L)
          ENDDO
          WQWDSL(L,KC,20) = XDSL(20) * CONV1
        ENDDO
      ENDIF
C
C *** C50 WET DEPOSTION (MULTIPLIED BY RAINFALL VOLUME IN RWQATM)
C
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      READ(1,90) TITLE(1)
      WRITE(2,999)
      DO M=1,3
        READ(1,999)
      ENDDO
      READ(1,*) (WQATM(NW),NW=1,6)
      WRITE(2,*) (WQATM(NW),NW=1,6)
      READ(1,*) (WQATM(NW),NW=7,13)
      WRITE(2,*) (WQATM(NW),NW=7,13)
      READ(1,*) (WQATM(NW),NW=14,NWQV)
      WRITE(2,*) (WQATM(NW),NW=14,NWQV)
      WRITE(2,999)
      WRITE(2,90) TITLE(1)
      WRITE(2, 21)' : (CHC, CHD, CHG)       = ',(WQATM(NW),NW=1,3)
      WRITE(2, 21)' : (ROC, LOC, DOC)       = ',(WQATM(NW),NW=4,6)
      WRITE(2, 21)' : (ROP, LOP, DOP, P4D)  = ',(WQATM(NW),NW=7,10)
      WRITE(2, 21)' : (RON, LON, DON)       = ',(WQATM(NW),NW=11,13)
      WRITE(2, 21)' : (NHX, NOX)            = ',(WQATM(NW),NW=14,15)
      WRITE(2, 21)' : (SUU, SAA, COD, DOX)  = ',(WQATM(NW),NW=16,19)
      WRITE(2,981)' : (TAM, FCB, CO2)       = ',(WQATM(NW),NW=20,NWQV)	!VB ADDED CO2
C
C *** C51 INPUT/OUTPUT FILE NAMES FOR SPATIALLY AND/OR TEMPORALLY VARYING PARAMETERS
C ***     MG/L FOR 1-19, TAM-MOLES/L, AND FCB-MPN/L
      IF(ISSKIP .GT. 0) CALL SKIPCOMM(1,CCMRM)
      READ(1,90) TITLE(1)
      WRITE(2,999)
      WRITE(2,90) TITLE(1)
      READ(1,295) RSTOFN
      WRITE(2,85)'* OUTPUT FILE FOR RESTART WRITING         = ', RSTOFN
      IF(IWQRST.EQ.1)THEN
        OPEN(99,FILE=RSTOFN,STATUS='UNKNOWN')
        CLOSE(99,STATUS='DELETE')
        OPEN(99,FILE=RSTOFN,STATUS='UNKNOWN')
        WRITE(99,888)
        CLOSE(99)
      ELSE
        IF(RSTOFN(1:4).NE.'NONE'.AND.RSTOFN(1:4).NE.'none')THEN
           PRINT*,ANS(PARTID2),RSTOFN(1:4)
           STOP 'ERROR!! INVALID IWQORST/RSTOFN'
         ENDIF
      ENDIF
  888 FORMAT('    L    K',
     &    '          BC          BD          BG        RPOC
     & LPOC','         DOC        RPOP        LPOP         DOP
     & PO4T','        RPON        LPON         DON         NH4
     & NO3','          SU          SA         COD          O2
     & TAM','         FCB			CO2			MAC')
      READ(1,295) ICIFN
      IF(IWQICI.EQ.1)THEN
        WRITE(2,"('SPATIALLY CONSTANT INITIAL CONDITIONS FROM C44')")
      ELSEIF(IWQICI.EQ.2)THEN
        WRITE(2,85)'* FILE FOR INITIAL CONDITIONS             = ', ICIFN
      ELSE
        IF(ICIFN(1:4).NE.'NONE'.AND.ICIFN(1:4).NE.'none')STOP 'ERROR!! INVALID IWQICI/ICIFN'
      ENDIF

      READ(1,295) AGRFN
      WRITE(2,85)'* FILE FOR ALGAL GROWTH, RESP., PREDATAT. = ', AGRFN
      IF(IWQAGR.NE.1)THEN
        IF(AGRFN(1:4).NE.'NONE'.AND.AGRFN(1:4).NE.'none')
     &     STOP 'ERROR!! INVALID IWQAGR/AGRFN'
      ENDIF

      READ(1,295) STLFN
      WRITE(2,85)'* FILE FOR SETTLING RATES OF ALGAE, PART. = ', STLFN
      IF(IWQSTL.NE.1)THEN
        IF(STLFN(1:4).NE.'NONE'.AND.STLFN(1:4).NE.'none')
     &     STOP 'ERROR!! INVALID IWQSTL/STLFN'
      ENDIF

      READ(1,295) SUNFN
      WRITE(2,85)'* FILE FOR IO, FD, TE, KT                 = ', SUNFN
      !IF(IWQSUN.NE.1)THEN
      !   IF(SUNFN(1:4).NE.'NONE'.AND.SUNFN(1:4).NE.'none')
      !&     STOP 'ERROR!! INVALID IWQSUN/SUNFN'
      !ENDIF

      READ(1,295) BENFN
      WRITE(2,85)'* FILE FOR BENTHIC FLUX                   = ', BENFN
      IF(IWQBEN.EQ.2)THEN
      ELSE
        IF(BENFN(1:4).NE.'NONE'.AND.BENFN(1:4).NE.'none')
     &     STOP 'ERROR!! INVALID IWQBEN/BENFN'
      ENDIF

      READ(1,295) PSLFN
      WRITE(2,85)'* FILE FOR POINT SOURCE INPUT             = ', PSLFN

      READ(1,295) NPLFN
      WRITE(2,85)'* FILE FOR NPS INPUT INCLUDING ATM. INPUT = ', NPLFN
      IF(IWQNPL.NE.1)THEN
        IF(NPLFN(1:4).NE.'NONE'.AND.NPLFN(1:4).NE.'none')
     &     STOP 'ERROR!! INVALID IWQNPL/NPLFN'
      ENDIF

      READ(1,295) NCOFN
      WRITE(2,85)'* DIAGNOSTIC FILE FOR NEGATIVE CONCENTRAT = ', NCOFN
      CLOSE(1)

      IF(IWQNC.EQ.1)THEN
        OPEN(1,FILE=NCOFN,STATUS='UNKNOWN')
        CLOSE(1,STATUS='DELETE')
        OPEN(1,FILE=NCOFN,STATUS='UNKNOWN')
      WRITE(1,284)'* NEGATIVE CONCENTRATION OCCURS:'
        CLOSE(1)
      ELSE
        IF(NCOFN(1:4).NE.'NONE'.AND.NCOFN(1:4).NE.'none')
     &     STOP 'ERROR!! INVALID IWQNC/NCOFN'
      ENDIF
  294 FORMAT(2I4,2I3, 7F8.3, /, 14X, 7F8.3, /, 14X, 9F8.3)
  295 FORMAT(44X, A50)
   96 FORMAT(2I5, 13I5, /, 10X, 8I5)
  969 FORMAT(2I4,1X,21I3)
  970 FORMAT(1X,21I3)
   97 FORMAT(2I4, 6F8.3, /, 8X, 7F8.3, /, 8X, 9F8.3)
   98 FORMAT(6F8.4, /, 7F8.4, /, 8F8.4)
   99 FORMAT(7F8.4, /, 7F8.4, /, 8F8.4)
   21 FORMAT(A27, 1P, 4E11.3)
  981 FORMAT(A27, 1P, 4E11.3)			!VB CHANGED 3E11.3 TO 4E11.3 TO INCLUDE CO2
   23 FORMAT(A46, I5)
   85 FORMAT(A44, A50)
  284 FORMAT(A32, /, 'NAME    ITNWQ    L    I    J    K       CONC')
      DO I=1,IWQZ
        IWQKA(I)=IWQKA(1)
        WQKRO(I)=WQKRO(1)
        WQKTR(I)=WQKTR(1)
        REAC(I)=REAC(1)
        WQKDC(I)=WQKDC(1)
        WQKDCALM(I)=WQKDCALM(1)
        WQKHRM(I)=WQKHRM(1)
        WQDOPM(I)=WQDOPM(1)
        WQKCD(I)=WQKCD(1)
        WQKHCOD(I)=WQKHCOD(1)
      ENDDO
      IF(IWQZ .GT. 1 .AND. IWQKIN .GT. 0)THEN
        PRINT *,'WQ: KINETICS.INP'
        OPEN(1,FILE='KINETICS.INP',STATUS='UNKNOWN')
        CALL SKIPCOMM(1,CCMRM)
        WRITE(2,*) ' '
        WRITE(2,*) ' SPATIALLY-VARYING KINETICS.INP FILE'
        WRITE(2,9111)
        DO I=1,IWQZ
          READ(1,*) IZ, IWQKA(IZ), WQKRO(IZ), WQKTR(IZ), REAC(IZ),
     &        WQKDC(IZ),WQKDCALM(IZ),WQKHRM(IZ),WQDOPM(IZ),WQKCD(IZ),
     &        WQKHCOD(IZ)
          WRITE(2,*) IZ, IWQKA(IZ), WQKRO(IZ), WQKTR(IZ), REAC(IZ),
     &        WQKDC(IZ),WQKDCALM(IZ),WQKHRM(IZ),WQDOPM(IZ),WQKCD(IZ),
     &        WQKHCOD(IZ)
          WRITE(2,9112) IZ, IWQKA(IZ), WQKRO(IZ), WQKTR(IZ), REAC(IZ),
     &        WQKDC(IZ),WQKDCALM(IZ),WQKHRM(IZ),WQDOPM(IZ),WQKCD(IZ),
     &        WQKHCOD(IZ)
        ENDDO
        CLOSE(1)
      ENDIF
 9111 FORMAT(/,'ZONE IWQKA   KRO   KTR  REAC   KDC KDCALGM  KHRM',
     &    ' DOPTM   KCD KHCOD')
 9112 FORMAT(I4, I6, 4F6.3, F8.3, 4F6.3)
C
C SET UP LOOK-UP TABLE FOR TEMPERATURE DEPENDENCY OVER WQTDMIN TO WQTDMAX
C
      DO M=1,NWQTD
C        WTEMP =1.00*REAL(M-1)*0.1 - 14.95
        WTEMP=WQTDTEMP(M)
        TT20 = WTEMP-20.0
        DO I=1,IWQZ
          WQKCOD(M,I) = WQKCD(I) * EXP( WQKTCOD*(WTEMP-WQTRCOD) )
          WQTDKR(M,I) = WQKTR(I)**TT20
          WRITE(2,2223)M,I,WQKTR(I),WQTDKR(M,I)
        ENDDO
      ENDDO
C
C READ IN MAPPING INFORMATION FOR SPATIALLY-VARYING PARAMETERS (UNIT #7)
C
      DO K=1,KC
        DO L=2,LA
          IWQZMAP(L,K)=1
        ENDDO
      ENDDO
      IF(IWQZ .GT. 1)THEN
        OPEN(1,FILE='WQWCMAP.INP',STATUS='UNKNOWN')
        WRITE(2,999)
        READ(1,30) (TITLE(M), M=1,3)
        WRITE(2,30) (TITLE(M), M=1,3)
C
C      READ(1,999)
C
        READ(1,999)
        WRITE(2,999)
        WRITE(2,32)
        IN=0
        IJKC=IC*JC*KC
        DO M=1,IJKC
          READ(1,*,END=1111) I,J,K,IWQZX
          I=XLOC(I)
          J=YLOC(J)
          L = LIJ(I,J)
          IN=IN+1
          IF (CELL_INSIDE_DOMAIN(L)) THEN
            IF(IJCT(I,J).LT.1 .OR. IJCT(I,J).GT.8)THEN
              PRINT*, 'I, J, K, IJCT(I,J) = ', I,J,K,IJCT(I,J)
              STOP 'ERROR!! INVALID (I,J) IN FILE WQWCMAP.INP'
            ENDIF
            IWQZMAP(L,K)=IWQZX
            WRITE(2,31) L,I,J,K,IWQZMAP(L,K)
          ENDIF
        ENDDO
 1111   CONTINUE
        IF(IN.NE.(LA-1)*KC)THEN
          PRINT*, 'ALL ACTIVE WATER CELLS SHOULD BE MAPPED FOR WQ PAR.'
      STOP 'ERROR!! NUMBER OF LINES IN FILE WQWCMAP.INP =\ (LA-1)'
      ENDIF
        CLOSE(1)
      ENDIF
C
C READ IN MAPPING INFORMATION FOR SPATIALLY-VARYING BENTHIC FLUXES.
C FORMULATED FOR PECONIC BAY DATA WHICH INCLUDES %MUD FOR EACH CELL AS
C WELL AS MAPPING TO BOTH MUD AND SAND FLUXES.  SUBROUTINE RWQBEN2
C CONTAINS THE CODE TO INTERPOLATE THE FINAL FLUX FOR THE CELL BASED
C ON THE PERCENT MUD AND THE MUD/SAND FLUXES.
C
      IF(IWQBEN .EQ. 2)THEN
        DO K=1,2
          DO L=2,LA
            IBENMAP(L,K)=1
            XBENMUD(L) = 0.50
          ENDDO
        ENDDO
        OPEN(1,FILE='WQBENMAP.INP',STATUS='UNKNOWN')
        WRITE(2,999)
        DO M=1,4
          READ(1,30) TITLE(M)
          WRITE(2,30) TITLE(M)
        ENDDO
C
C SKIP ALL COMMENT CARDS AT BEGINNING OF FILE:
C
        REWIND(1)
        CCMRM = '#'
        CALL SKIPCOMM(1, CCMRM)
C
C       READ(1,999)
C
        WRITE(2,999)
        WRITE(2,33)
        IN=0
        IJKC=IC*JC
        DO M=1,IJKC
          READ(1,*,END=1112) I, J, XMUD, IZMUD, IZSAND
          I=XLOC(I)
          J=YLOC(J)
          L = LIJ(I,J)
          IN=IN+1
          IF (CELL_INSIDE_DOMAIN(L)) THEN
            IF(IJCT(I,J).LT.1 .OR. IJCT(I,J).GT.8)THEN
              PRINT*, 'I, J, K, IJCT(I,J) = ', I,J,IJCT(I,J)
              STOP 'ERROR!! INVALID (I,J) IN FILE WQBENMAP.INP'
            ENDIF
            IBENMAP(L,1) = IZMUD
            IBENMAP(L,2) = IZSAND
            XBENMUD(L) = XMUD / 100.0
            WRITE(2,34) L, I, J, XBENMUD(L), IBENMAP(L,1), IBENMAP(L,2)
          ENDIF
        ENDDO
 1112   CONTINUE
        IF(IN .NE. (LA-1))THEN
          PRINT*, 'ALL ACTIVE WATER CELLS SHOULD BE MAPPED FOR WQ PAR.'
          STOP 'ERROR!! NUMBER OF LINES IN FILE WQBENMAP.INP <> (LA-1)'
        ENDIF
        CLOSE(1)
      ENDIF
      CLOSE(2)
 2222 FORMAT(' M,WQKTR(1),WQTDKR(M,1) = ',I5,2F10.4)
 2223 FORMAT(' M,I,WQKTR(1),WQTDKR(M,I) = ',2I5,2F10.4)
   30 FORMAT(A79)
   31 FORMAT(15I5)
   32 FORMAT('    L    I    J    L IWQZMAP')
   33 FORMAT('     L    I    J   MUD IZMUD IZSAND')
   34 FORMAT(' ',3I5, F6.2, 2I6)
      RETURN
      END
